* ======================================================================
*
        PROGRAM NONH
*
*                C O P Y R I G H T -- 1994
*
*     A. HIBBERT, DEP'T OF APPLIED MATHEMATICS AND THEORETICAL PHYSICS
*        THE QUEENS UNIVERISTY OF BELFAST
*
*     C. FROESE FISCHER, DEP'T OF COMPUTER SCIENCE
*        VANDERBILT UNIVERISTY
*
*     Spectroscopic Notation developed by:
*     K. M. S. SAXENA, DEP'T OF THEORETICAL CHEMISTRY
*        UNIVERSITY OF ALBERTA
*
*     AUGUST, 1982
*
*     Computer Physics Communications, Vol. 64, 417--430 (1991).
* ======================================================================
*     GENERAL PROGRAM TO COMPUTE MATRIX ELEMENTS OF THE HAMILTONIAN.
*     UNDER THE FOLLOWING ASSUMPTIONS -
*         (1) LS-COUPLING
*         (2) ORTHO-NORMAL CONFIGURATION STATE FUNCTIONS
*         (3) AT MOST TWO (2) NON-ORTHOGONAL ORBITAL PAIRS
*             BETWEEN THE LEFT AND RIGHT CONFIGURATION STATES
*         (4) ALL ORBITALS WITHIN A PAIR ARE ORTHOGONAL
*
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      PARAMETER (NWD=30,NCD=100,NORD=NWD*(NWD-1)/2)
*
      COMMON/DEBUG/IBUG1,IBUG2,IBUG3,NBUG6,NBUG7,IFULL
      COMMON/INFORM/IREAD,IWRITE,IOUT,ISC0,ISC1,ISC2,ISC3,JSC(3),IALL
      COMMON/STATES/NCFG,MAXORB,IAJCMP(NWD),LJCOMP(NWD),
     :NJCOMP(NWD),NOCCSH(NCD),NELCSH(5,NCD),NOCORB(5,NCD),J1QNRD(9,NCD)
      COMMON /DIAGNL/IDIAG,JA,JB
      COMMON /FOUT/NOV(2),IOVLAP(10,2),NF,NG,NR,NL,NZ,NN,NV,NS,IFLAG,NIJ
          CHARACTER*2 ANS, INPUT*24
*
    1 FORMAT(//' IOUT =  INT.LST (OUTPUT FILE)'//)
*
    2 FORMAT(///20X,'============================'/
     :          20X,'   N O N - W E I G H T S    '/
     :          20X,'============================'/)
    3 FORMAT(/ ' IBUG1  =',I3,' (DEBUG IN WEIGHTS OF 1-EL PART)'/
     :         ' IBUG2  =',I3,' (DEBUG IN WEIGHTS OF 2-EL PART)'/
     :         ' IBUG3  =',I3,' (DEBUG IN RECOUPLING PACKAGE)'//)
*
*     ...  THE FOLLOWING SECTION CONCERNS INPUT/OUTPUT AND MAY BE
*          SYSTEM DEPENDENT.  CHECK ALLOWED UNIT NUMBERS AND
*          FILE NAME CONVENTIONS - MODIFY, IF NECESSARY.
*
CSUN      i = iargc()
CSUN      if (i .eq. 0) then
                 INPUT = 'cfg.inp'
CSUN      else
CSUN             call getarg(1,INPUT)
CSUN      end if
      IREAD=4
      IWRITE = 6
      IOUT=7
      ISC0=10
      ISC1=11
      ISC2=12
      ISC3=13
      WRITE(IWRITE,2)
*
      OPEN(UNIT=IOUT, FILE='int.lst',STATUS='UNKNOWN')
      OPEN(UNIT=ISC0,STATUS='SCRATCH',FORM='UNFORMATTED')
      OPEN(UNIT=ISC1,STATUS='SCRATCH',FORM='UNFORMATTED')
      OPEN(UNIT=ISC2,STATUS='SCRATCH',FORM='UNFORMATTED')
      OPEN(UNIT=ISC3,STATUS='SCRATCH',FORM='UNFORMATTED')
*     OPEN(UNIT=ISC0,FILE='isc0',STATUS='SCRATCH')
*     OPEN(UNIT=ISC1,FILE='isc1',STATUS='SCRATCH')
*     OPEN(UNIT=ISC2,FILE='isc2',STATUS='SCRATCH')
*     OPEN(UNIT=ISC3,FILE='isc3',STATUS='SCRATCH')
*
*     ... END OF MACHINE DEPENDENT SECTION
*
      IBUG1 = 0
      IBUG2 = 0
      IBUG3 = 0
*
*     ... The following statements allow the setting of DEBUG
*
CDBG  WRITE(0,'(A)') ' IBUG1,IBUG2,IBUG3 ? '
CDBG  READ(5,*) IBUG1,IBUG2,IBUG3
CDBG  WRITE(IWRITE,1) IBUG1,IBUG2,IBUG3
      WRITE(0,'(A)') ' FULL PRINT-OUT ? (Y/N) '
      IFULL = 0
      READ(5,'(A2)') ANS
      IF ( ANS .EQ. 'Y' .OR. ANS .EQ. 'y') IFULL = 1
*
      NEW = 0
      NZER0 = 0
      NONE = 0
      IALL = 1
*
*  ---  Note: IALL=0 allows the output of an INT.LST to be in the
*   format required by MCHF88, but SUBROUTINE OUTLS must then be
*   simplified to avoid the sort and output the Fk, Gk, Rk and L
*   data in the old format.
*
*  ---  Determine input data; non-orthogonal case
*
      CALL CFGN1(INPUT)
*
      WRITE(0,'(/A)') ' ALL INTERACTIONS ? (Y/N) '
      READ(5,'(A1)') ANS
      IF (ANS .EQ. 'N' .OR. ANS .EQ. 'n' ) THEN
          WRITE(0,'(A)') ' NEW (0=ALL), N-ZERO (0=ALL) ? '
          READ(5,*) NEW, NZERO
      END IF
      IF ( NEW .EQ. 0 ) NEW = NCFG
          IF (NZERO .EQ. 0) THEN
                 NZERO=NCFG
                 IFIRST=0
          ELSE
                 IFIRST=1
          END IF
*
*     SET FACTORIALS AND LOG OF FACTORIALS
*
      CALL INITA
*
*     Initialize parameters for output
      NOV(1) = 0
      NOV(2) = 0
      NF = 0
      NG = 0
      NR = 0
      NL = 0
      NIJ = 0
*
      CALL ANGMOM(NEW,NZERO,IFIRST)
      CALL OUTLS
 10   CONTINUE
      STOP
      END
*
*     ------------------------------------------------------------------
*       A L L A D D
*     ------------------------------------------------------------------
*
      SUBROUTINE ALLADD(IHSH,M3,M4,M5,M6,M7,M8,M9,M10,
     :M11,M12,M13,M14,M15,M16,M17,M18)
*
* --- Sets up variables related to IHSH
*
      M3=IHSH-1
      M4=IHSH+1
      M5=IHSH+2
      M6=2*IHSH-1
      M7=M6+1
      M8=M3+M6
      M9=M8+1
      M10=M8+2
      M11=M8+3
      M12=M8+4
      M13=M8+5
      M14=M8+6
      M15=M8+7
      M16=M8+8
      M17=M8+9
      M18=IHSH+3
      RETURN
      END
*
*     ------------------------------------------------------------------
*       A N G M O M
*     ------------------------------------------------------------------
*
      SUBROUTINE ANGMOM(NEW,NZERO,IFIRST)
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      PARAMETER (NWD=30,NCD=100,NORD=NWD*(NWD-1)/2)
*
      COMMON/INFORM/IREAD,IWRITE,IOUT,ISC0,ISC1,ISC2,ISC3,JSC(3),IALL
      COMMON/DEBUG/IBUG1,IBUG2,IBUG3,NBUG6,NBUG7,IFULL
      COMMON/OVRLAP/MU,NU,MUP,NUP,NONORT,NOVLPS,IROWMU,IROWNU,ICOLMU,
     : ICOLNU,NORTH,IORDER,NCALLS,LMU,LNU,LMUP,LNUP,JMU,JNU,JMUP,JNUP,
     :     IORTH(NORD)
      COMMON/STATES/NCFG,MAXORB,IAJCMP(NWD),LJCOMP(NWD),
     :NJCOMP(NWD),NOCCSH(NCD),NELCSH(5,NCD),NOCORB(5,NCD),J1QNRD(9,NCD)
      COMMON/DIAGNL/IDIAG,JA,JB
      COMMON /FOUT/NOV(2),IOVLAP(10,2),NF,NG,NR,NL,NZ,NN,NV,NS,IFLAG,NIJ
      CHARACTER*30 FORMAT(3)
      CHARACTER*1 NCHAR(5)
      DATA NCHAR/'1','2','3','4','5'/
      DATA FORMAT/'(2H <,I3,5H |H| ,I2,6H > = <, ',
     :              '  (A3,1H(,I2,1H)) ,5H |H| ,   ',
     :                '(A3,1H(,I2,1H)),2H >,/)     '/
*
* --- THIS PROGRAMME CONSIDERS EITHER SUPERPOSITION OF CONFIGURATIONS OR
*     MULTI-CONFIGURATIONAL HARTREE-FOCK WAVE FUNCTIONS.  USING THE
*     RESULT THAT THE TWO-ELECTRON HAMILTONIAN MATRIX ELEMENT
*     (PSI/V/PSIP)  CAN BE WRITTEN AS A SUM OF SLATER INTEGRALS, THE
*     PRESENT CODE  -  WEIGHTS  -  CALCULATES THE COEFFICIENTS OF THESE
*     INTEGRALS.  PSI AND PSIP ARE ALLOWED TO RUN OVER NCFG CONFIGURATNS
*
*
* --- CONSIDER (PSI/V/PSIP) AS PSI AND PSIP RUN OVER ALL CONFIGURATIONS
*
      NFIRST = NCFG - NEW + 1
      DO 6 JA=NFIRST,NCFG
      WRITE(0,*) '   ',ja
      DO 6 JB=1,JA
      IF (JB .GT. NZERO .AND. IFIRST .EQ. 1 .AND. JA .NE. JB ) GO TO 6
      IF (JB.GT.NZERO .AND. JB.LT.NFIRST .AND. IFIRST.EQ.0) GO TO 6
      IFLAG = 0
      IDIAG=0
      IF(JA.EQ.JB) IDIAG=1
      IF (NORTH .NE. 0) CALL NORTBP(JA,JB)
      N1=NOCCSH(JA)
      N2=NOCCSH(JB)
      IF (IFULL .NE. 0) THEN
         FORMAT(2)(2:2) = NCHAR(N1)
         FORMAT(2)(30:30) = NCHAR(N2)
         WRITE(IWRITE,'(///)')
         WRITE(IWRITE,FORMAT) JA,JB,
     :        (IAJCMP(NOCORB(J,JA)),NELCSH(J,JA),J=1,N1),
     :        (IAJCMP(NOCORB(J,JB)),NELCSH(J,JB),J=1,N2)
      END IF
*
* --- SET UP DEFINING QUANTUM NUMBERS FOR EACH MATRIX ELEMENT
*
      CALL SETUP(JA,JB)
      IF(IBUG1.GT.0.OR.IBUG2.GT.0) CALL VIJOUT(JA,JB)
*
* --- TEST ON POSSIBLE RECOUPLING ORTHOGONALITY
*
      CALL ORTHOG(LET)
      IF(LET.EQ.0) GO TO 6
*
* --- IF NO SUCH ORTHOGONALITY IS EXHIBITED, CALCULATE WEIGHTS OF SLATER
*     INTEGRALS
*
      CALL H0WTS
      CALL CHOP
      CALL RKWTS
      IF (IFLAG .NE. 0) NIJ = NIJ + 1
    6 CONTINUE
      NTOTAL = ((2*NCFG - NEW +1)*NEW)/2
      WRITE(0,*) NTOTAL, ' matrix elements'
      WRITE(0,*) NIJ, ' non-zero matrix elements'
      WRITE(0,*) 100*NIJ/REAL(NTOTAL),' % dense'
      WRITE(0,*) 'NF=',nf,' NG=',ng,' NR=',nr,' NL=',nl
      WRITE(0,*) 'Total number of terms =', NF+NG+NR+NL
 1011 RETURN
      END
*
*     ------------------------------------------------------------------
*       C H O P
*     ------------------------------------------------------------------
*
      SUBROUTINE CHOP
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      COMMON/MEDEFN/IHSH,NJ(10),LJ(10),NOSH1(10),NOSH2(10),J1QN1(19,3),
     :     J1QN2(19,3),IJFUL(10)
      COMMON/REMOVE/ICHOP(10)
      COMMON/DIAGNL/IDIAG,JA,JB
*
* --- ZEROIZE THE OUTPUT ARRAY
*
      DO 1 I=1,IHSH
      ICHOP(I)=0
    1 CONTINUE
*
*     NO AVERAGE ENERGY TERMS FOR OFF-DIAGONAL MATRIX ELEMENTS
*
      IF(IDIAG.EQ.0) RETURN
      JSTO=0
      ICOUNT=0
      DO 3 J=1,IHSH
      NFULL=4*LJ(J)+2
      I2=NOSH1(J)
*
*     IS THE SHELL FULL OR EMPTY
*
      IF(I2.EQ.NFULL.OR.I2.EQ.0) GO TO 4
*
*     IF NOT, DOES IT CONTAIN ONLY ONE ELECTRON, OR ONLY ONE =HOLE=
*
      IF(I2.EQ.1.OR.I2.EQ.(NFULL-1)) JSTO=J
      GO TO 3
    4 ICOUNT=ICOUNT+1
*
*     ICHOP  SET UNITY FOR CLOSED SHELLS
*
      ICHOP(J)=1
    3 CONTINUE
*
*     IF ALL BUT ONE SHELL IS CLOSED, AND THIS CONTAINS ONE ELECTRON OR
*     =HOLE= , THEN IT CAN BE TREATED PURELY BY AVERAGE ENERGY
*
      IF(ICOUNT.NE.(IHSH-1).OR.JSTO.EQ.0) RETURN
      ICHOP(JSTO)=1
      RETURN
      END
*
*     ------------------------------------------------------------------
*       C L S H E L
*     ------------------------------------------------------------------
*
      SUBROUTINE CLSHEL (ISIG, ISIGP, IRHO )
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      PARAMETER (NWD=30,NCD=100,NORD=NWD*(NWD-1)/2)
*
      COMMON /CLOSED/B1ELC(4),NCLOSD,IAJCLD(NWD),LJCLSD(NWD),IBK
      COMMON/ENAV/NINTS,KVALUE(15),COEFCT(15)
      COMMON/STATES/NCFG,MAXORB,IAJCMP(NWD),LJCOMP(NWD),
     :NJCOMP(NWD),NOCCSH(NCD),NELCSH(5,NCD),NOCORB(5,NCD),J1QNRD(9,NCD)
      COMMON/MEDEFN/IHSH,NJ(10),LJ(10),NOSH1(10),NOSH2(10),J1QN1(19,3),
     :     J1QN2(19,3),IJFUL(10)
      COMMON/INFORM/IREAD,IWRITE,IOUT,ISC0,ISC1,ISC2,ISC3,JSC(3),IALL
      COMMON/DEBUG/IBUG1,IBUG2,IBUG3,NBUG6,NBUG7,IFULL
      COMMON/DIAGNL/IDIAG,JA,JB
      COMMON/SIGNF/SIGNFA
      COMMON/OVRINT/IOVEL(2),IOVER(2),IOVEP(2)
      COMMON/OVRLAP/MU,NU,MUP,NUP,NONORT,NOVLPS,IROWMU,IROWNU,ICOLMU,
     : ICOLNU,NORTH,IORDER,NCALLS,LMU,LNU,LMUP,LNUP,JMU,JNU,JMUP,JNUP,
     :     IORTH(NORD)
      DIMENSION IOVL(2),IOVR(2)
*
* --- Evaluates Hamiltonian matrix element contributions where one 
*     or more interacting shells are closed
*
  101 FORMAT(3F14.6,4X,1HF,I2,1H(,A3,1H,,A3,1H))
  102 FORMAT(3F14.6,4X,1HG,I2,1H(,A3,1H,,A3,1H))
  103 FORMAT(3F14.6,4X,1HR,I2,1H(,A3,1H,,A3,1H;,A3,1H,,A3,1H),
     :       :2(2H <,A3,1H|,A3,3H>**,I2:))
  201 FORMAT(F14.8,1HF,I2,1H(,A3,I3,1H,,A3,I3,1H))
  202 FORMAT(F14.8,1HG,I2,1H(,A3,I3,1H,,A3,I3,1H))
  203 FORMAT(F14.8,1HR,I2,1H(,2A3,I3,1H,,2A3,I3,1H),
     :      :2(1H<,A3,1H|,A3,1H>,I2:))
      ZERO = 0.D0
      IZERO = 0
      IBK1 = IBK+1
      IBK2 = IBK+2
      IF (IRHO .NE. 0) GO TO 20
      IF (ISIG .NE. 0) GO TO 1
*
* --- DIAGONAL CASE - ADD CLOSED SHELL CONTRIBUTIONS
*
      DO 2 J = 1,NCLOSD
      LA = LJCLSD(J)
      NA = 4*LA + 2
*
* --- CLOSED-CLOSED INTERACTIONS
*
      DO 3 K = 1,J
      LB = LJCLSD(K)
      NB = 4*LB + 2
      IF (J .EQ. K) GO TO 4
      IEQUIV = 2
      AC2 = DFLOAT(NA*NB)
      GO TO 5
    4 IEQUIV = 1
      AC2 = DFLOAT((NA*(NA-1))/2)
    5 CALL INTACT(LA,LB,IEQUIV)
      IF (IFULL .NE. 0)
     :   WRITE(IWRITE,101) AC2,AC2,ZERO,IZERO,IAJCLD(K),IAJCLD(J)
      IF(NINTS .EQ. 0) GO TO 3
      DO 7 N = 1,NINTS
      ZA = AC2*COEFCT(N)
      K1 = KVALUE(N)
      IF (IEQUIV .EQ. 1) GO TO 8
      IF (IFULL .NE. 0)
     :   WRITE (IWRITE, 102) ZA,ZA, ZERO,K1,IAJCLD(K),IAJCLD(J)
      GO TO 7
    8 IF (IFULL .NE. 0)
     :   WRITE(IWRITE, 101) ZA,ZA,ZERO,K1,IAJCLD(K),IAJCLD(J)
    7 CONTINUE
    3 CONTINUE
*
* --- CLOSED - OPEN INTERACTIONS
*
      DO 9 K = 1, IHSH
      LB = LJ(K)
      NB = NOSH1(K)
      JFULL = IJFUL(K)
      IEQUIV = 2
      AC2 = DFLOAT(NA*NB)
      CALL INTACT(LA,LB,IEQUIV)
      IF (IFULL .NE. 0)
     :   WRITE (IWRITE, 101) AC2,AC2,ZERO,IZERO,IAJCLD(J),IAJCMP(JFULL)
      IF (NINTS .EQ. 0) GO TO 9
      DO 10 N = 1,NINTS
      ZA = AC2*COEFCT(N)
      K1 = KVALUE(N)
      IF (IFULL .NE. 0)
     :   WRITE (IWRITE, 102) ZA, ZA,ZERO,K1,IAJCLD(J),IAJCMP(JFULL)
   10 CONTINUE
    9 CONTINUE
    2 CONTINUE
      RETURN
*
* --- OFF-DIAGONAL, ONE-ELECTRON DIFFERENT, ADD COMMON
*     CLOSED-SHELL CONTRIBUTIONS
*
    1 IF (DABS(B1ELC(IBK1))+DABS(B1ELC(IBK2)) .LT. 1.D-10) RETURN
*
* --- SET UP THE OVERLAPS
*
      IF (NOVLPS .GE. 1) THEN
         IOVL(1) = IJFUL(IOVEL(1))
         IOVR(1) = IJFUL(IOVER(1))
         IF (NOVLPS .EQ. 2) THEN
            IOVL(2) = IJFUL(IOVEL(2))
            IOVR(2) = IJFUL(IOVER(2))
         ENDIF
       ENDIF
      LB = LJ(ISIG)
      JSIG = IJFUL(ISIG)
      JSIGP = IJFUL(ISIGP)
      IEQUIV = 2
      DO 13 IOV = 1,NONORT+1
      IF (DABS(B1ELC(IOV+IBK)) .LT. 1.D-10) GO TO 14
      DO 11 J = 1,NCLOSD
      LA = LJCLSD(J)
      NA = 4*LA + 2
      AC2 = B1ELC(IOV+IBK)*NA*SIGNFA
      CALL INTACT(LA,LB,IEQUIV)
      IF (IFULL .NE. 0)
     :    WRITE(IWRITE,103)AC2,ZERO,AC2,IZERO,IAJCLD(J),IAJCMP(JSIG),
     :   IAJCLD(J),IAJCMP(JSIGP),
     :   (IAJCMP(IOVL(II)),IAJCMP(IOVR(II)),IOVEP(II),II=1,NOVLPS)
      IF (IALL.EQ.0) THEN
         CALL OVLSET(NOVLPS,IOVL(1),IOVR(1),IOVEP(1),
     :                      IOVL(2),IOVR(2),IOVEP(2),IPTR)
         CALL SAVE(3,2.*AC2,IZERO,64-J,JSIG,64-J,JSIGP,JA,JB,IPTR)
      END IF
   15 IF (NINTS .EQ. 0) GO TO 11
      DO 12 N=1,NINTS
      ZA = AC2*COEFCT(N)
      K = KVALUE(N)
      IF (IFULL .NE. 0) WRITE(IWRITE,103) ZA,ZERO,ZA,K,IAJCLD(J),
     :    IAJCMP(JSIG),IAJCMP(JSIGP),IAJCLD(J),
     :    (IAJCMP(IOVL(II)),IAJCMP(IOVR(II)),IOVEP(II),II=1,NOVLPS)
      IF (IALL.EQ.0) THEN
         CALL OVLSET(NOVLPS,IOVL(1),IOVR(1),IOVEP(1),
     :                      IOVL(2),IOVR(2),IOVEP(2),IPTR)
         CALL SAVE(3,2.*ZA,K,64-J,JSIG,JSIGP,64-J,JA,JB,IPTR)
      END IF
   12 CONTINUE
   11 CONTINUE
   14 IOTEMP = IOVR(1)
      IOVR(1) = IOVR(2)
      IOVR(2) = IOTEMP
   13 CONTINUE
      RETURN
*
* --- OFF-DIAGONAL, ONE-ELECTRON DIFFERENT, TREAT BY AVERAGE
*     ENERGY FORMULAE
*
   20 IF (DABS(B1ELC(IBK+1))+DABS(B1ELC(IBK+2)) .LT. 1.D-10) RETURN
*
* --- SET UP THE OVERLAPS
*
      IF (NOVLPS .GE. 1) THEN
         IOVL(1) = IJFUL(IOVEL(1))
         IOVR(1) = IJFUL(IOVER(1))
         IF (NOVLPS .EQ. 2) THEN
            IOVL(2) = IJFUL(IOVEL(2))
            IOVR(2) = IJFUL(IOVER(2))
         ENDIF
      ENDIF
      LB = LJ(ISIG)
      LA = LJ(IRHO)
      JSIG = IJFUL(ISIG)
      JSIGP = IJFUL(ISIGP)
      JRHO = IJFUL(IRHO)
      IEQUIV = 2
      CALL INTACT(LA,LB,IEQUIV)
      NA = 4*LA + 2
      DO 23 IOV = 1,NONORT+1
      AC2 = B1ELC(IOV+IBK)*NA*SIGNFA
      IF (IFULL .NE. 0)
     :   WRITE(IWRITE,103) AC2,ZERO,AC2,IZERO,IAJCMP(JRHO),IAJCMP(JSIG),
     :   IAJCMP(JRHO),IAJCMP(JSIGP),
     :   (IAJCMP(IOVL(II)),IAJCMP(IOVR(II)),IOVEP(II),II=1,NOVLPS)
      Y = AC2
      IF (IALL .EQ. 0) Y = 2.D0*Y
      CALL OVLSET(NOVLPS,IOVL(1),IOVR(1),IOVEP(1),
     :                   IOVL(2),IOVR(2),IOVEP(2),IPTR)
      CALL SAVE(3,Y,IZERO,JRHO,JSIG,JRHO,JSIGP,JA,JB,IPTR)
   22 IF (NINTS .EQ. 0) GO TO 24
      DO 21 N=1,NINTS
      ZA = AC2*COEFCT(N)
      K = KVALUE(N)
      IF (IFULL .NE. 0)
     :    WRITE (IWRITE,103) ZA,ZERO,ZA,K,IAJCMP(JRHO),IAJCMP(JSIG),
     :    IAJCMP(JSIGP),IAJCMP(JRHO),
     :    (IAJCMP(IOVL(II)),IAJCMP(IOVR(II)),IOVEP(II),II=1,NOVLPS)
      Y = ZA
      IF (IALL .EQ. 0) Y = 2.D0*Y
      CALL OVLSET(NOVLPS,IOVL(1),IOVR(1),IOVEP(1),
     :                   IOVL(2),IOVR(2),IOVEP(2),IPTR)
      CALL SAVE(3,Y,K,JRHO,JSIG,JSIGP,JRHO,JA,JB,IPTR)
   21 CONTINUE
   24 IOTEMP = IOVR(1)
      IOVR(1) = IOVR(2)
      IOVR(2) = IOTEMP
   23 CONTINUE
      RETURN
      END
*
*     ------------------------------------------------------------------
*       C N D E N S
*     ------------------------------------------------------------------
*
      SUBROUTINE CNDENS(I1L)
*
      PARAMETER (NWD=30,NCD=100,NORD=NWD*(NWD-1)/2)
      PARAMETER(KFL1=60,KFL2=12)
      LOGICAL FREE
      COMMON/COUPLE/NJ1S,NJ23S,J1(KFL1),J2(KFL2,3),J3(KFL2,3),FREE(KFL1)
      COMMON/OVRLAP/MU,NU,MUP,NUP,NONORT,NOVLPS,IROWMU,IROWNU,ICOLMU,
     : ICOLNU,NORTH,IORDER,NCALLS,LMU,LNU,LMUP,LNUP,JMU,JNU,JMUP,JNUP,
     :     IORTH(NORD)
*
* --- Links the co-ordinates of electrons in non-orthogonal subshells 
*     on the two sides of the martix element, and condenses the angular 
*     momentum coupling schemes accordingly.
*
      NJ23ST=NJ23S
      NJ23=NJ23S-1
      CALL PIKUP1(J2,NJ23,MUP,NUP,IROWMU,IROWNU,ICOLMU,ICOLNU)
      IF(IROWMU.NE.IROWNU) GO TO 1
      IROWM1=IROWMU+1
      J2(IROWM1+1,1)=J2(IROWM1,2)
      CALL SCRAP(J2,NJ23,IROWMU)
      CALL SCRAP(J2,NJ23,IROWMU)
      GO TO 3
    1 IMU=J2(IROWMU,3)
      CALL PIKUP2(J2,NJ23,IROWMU,IMU,IROW,ICOL)
      IF(IROWMU.EQ.NJ23) GO TO 2
      J2(IROW,ICOL)=J2(IROWMU,3-ICOLMU)
    2 IF(NOVLPS.EQ.1.OR.MUP.EQ.NUP) GO TO 4
      IMU=J2(IROWNU,3)
      CALL PIKUP2(J2,NJ23,IROWNU,IMU,IROW,ICOL)
      IF(IROWNU.EQ.NJ23) GO TO 4
      J2(IROW,ICOL)=J2(IROWNU,3-ICOLNU)
    4 CALL SCRAP(J2,NJ23,IROWMU)
      IF(NOVLPS.LT.2.OR.MUP.EQ.NUP) GO TO 3
      IROWNU=IROWNU-1
      CALL SCRAP(J2,NJ23,IROWNU)
    3 NJ23=NJ23ST-1
      CALL PIKUP1(J3,NJ23,MU,NU,IROWMU,IROWNU,ICOLMU,ICOLNU)
      IF(IROWMU.NE.IROWNU) GO TO 11
      IROWM1=IROWMU+1
      J3(IROWM1+1,1)=J3(IROWM1,2)
      CALL SCRAP(J3,NJ23,IROWMU)
      CALL SCRAP(J3,NJ23,IROWMU)
      GO TO 13
   11 IMU=J3(IROWMU,3)
      CALL PIKUP2(J3,NJ23,IROWMU,IMU,IROW,ICOL)
      IF(IROWMU.EQ.NJ23) GO TO 12
      J3(IROW,ICOL)=J3(IROWMU,3-ICOLMU)
   12 IF(NOVLPS.EQ.1.OR.MU.EQ.NU) GO TO 14
      IMU=J3(IROWNU,3)
      CALL PIKUP2(J3,NJ23,IROWNU,IMU,IROW,ICOL)
      IF(IROWNU.EQ.NJ23) GO TO 14
      J3(IROW,ICOL)=J3(IROWNU,3-ICOLNU)
   14 CALL SCRAP(J3,NJ23,IROWMU)
      IF(NOVLPS.LT.2.OR.MU.EQ.NU) GO TO 13
      IROWNU=IROWNU-1
      CALL SCRAP(J3,NJ23,IROWNU)
   13 IF(MU.EQ.NU.AND.MUP.NE.NUP) GO TO 21
      IF(MUP.EQ.NUP.AND.MU.NE.NU) GO TO 22
      CALL PIKUP1(J3,NJ23,MUP,NUP,IROWMU,IROWNU,ICOLMU,ICOLNU)
      J3(IROWMU,ICOLMU)=MU
      IF(IROWNU.EQ.0) GO TO 23
      J3(IROWNU,ICOLNU)=NU
      GO TO 23
   21 NJ23=NJ23ST-3
      IROW=3-I1L
      CALL NSCRAP(J2,NJ23,IROW)
      J2(IROW,1)=MUP
      J2(IROW,2)=NUP
      J2(IROW,3)=MU
      GO TO 23
   22 NJ23=NJ23ST-3
      IROW=3-I1L
      CALL NSCRAP(J3,NJ23,IROW)
      J3(IROW,1)=MU
      J3(IROW,2)=NU
      J3(IROW,3)=MUP
   23 NJ23S=NJ23+1
      RETURN
      END
*
*     ------------------------------------------------------------------
*       F A N O
*     ------------------------------------------------------------------
*
      SUBROUTINE FANO(IRHO,ISIG,IRHOP,ISIGP)
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      PARAMETER (NWD=30,NCD=100,NORD=NWD*(NWD-1)/2)
      DIMENSION RMEDIR(19),RMEEX(19),NBAR(10)
 
      PARAMETER(KFL1=60,KFL2=12,
     +          KFL6=120,KFL7=150,KFL8=120,KFL9=40,KFLW=20,
     +          KFLS=12,KFLN=10,KFLV=40)
C
      LOGICAL FAILSD1,FAILSE1,FAILAD1,FAILAE1,FREE
C
      LOGICAL SMVRSD1,SMVRSE1,SMVRAD1,SMVRAE1
      DIMENSION K6SD1(KFL6),K7SD1(KFL7),K8SD1(KFL8),K9SD1(KFL9),
     :          KWSD1(6,KFLW),LDELSD1(KFLW,2),SMVRSD1(KFL1)
      DIMENSION K6SE1(KFL6),K7SE1(KFL7),K8SE1(KFL8),K9SE1(KFL9),
     :          KWSE1(6,KFLW),LDELSE1(KFLW,2),SMVRSE1(KFL1)
      DIMENSION K6AD1(KFL6),K7AD1(KFL7),K8AD1(KFL8),K9AD1(KFL9),
     :          KWAD1(6,KFLW),LDELAD1(KFLW,2),SMVRAD1(KFL1)
      DIMENSION K6AE1(KFL6),K7AE1(KFL7),K8AE1(KFL8),K9AE1(KFL9),
     :          KWAE1(6,KFLW),LDELAE1(KFLW,2),SMVRAE1(KFL1)
C
      DIMENSION J6PSD1(KFLV),J7PSD1(KFLV),J8PSD1(KFLV),J9PSD1(KFLV),
     :    JWRDSD1(6,KFLW),NBJSD1(KFLN),NB6JSD1(KFLN),K6CPSD1(KFLN),
     :    K7CPSD1(KFLN),K8CPSD1(KFLN),K9CPSD1(KFLN),JSM6SD1(KFLS),
     :    JSM4SD1(KFLS,KFLW),JSM5SD1(KFLS,KFLW),IN6JSD1(KFLW)
      DIMENSION J6PSE1(KFLV),J7PSE1(KFLV),J8PSE1(KFLV),J9PSE1(KFLV),
     :    JWRDSE1(6,KFLW),NBJSE1(KFLN),NB6JSE1(KFLN),K6CPSE1(KFLN),
     :    K7CPSE1(KFLN),K8CPSE1(KFLN),K9CPSE1(KFLN),JSM6SE1(KFLS),
     :    JSM4SE1(KFLS,KFLW),JSM5SE1(KFLS,KFLW),IN6JSE1(KFLW)
      DIMENSION J6PAD1(KFLV),J7PAD1(KFLV),J8PAD1(KFLV),J9PAD1(KFLV),
     :    JWRDAD1(6,KFLW),NBJAD1(KFLN),NB6JAD1(KFLN),K6CPAD1(KFLN),
     :    K7CPAD1(KFLN),K8CPAD1(KFLN),K9CPAD1(KFLN),JSM6AD1(KFLS),
     :    JSM4AD1(KFLS,KFLW),JSM5AD1(KFLS,KFLW),IN6JAD1(KFLW)
      DIMENSION J6PAE1(KFLV),J7PAE1(KFLV),J8PAE1(KFLV),J9PAE1(KFLV),
     :    JWRDAE1(6,KFLW),NBJAE1(KFLN),NB6JAE1(KFLN),K6CPAE1(KFLN),
     :    K7CPAE1(KFLN),K8CPAE1(KFLN),K9CPAE1(KFLN),JSM6AE1(KFLS),
     :    JSM4AE1(KFLS,KFLW),JSM5AE1(KFLS,KFLW),IN6JAE1(KFLW)
C
      LOGICAL FAILSD2,FAILSE2,FAILAD2,FAILAE2
C
      LOGICAL SMVRSD2,SMVRSE2,SMVRAD2,SMVRAE2
      DIMENSION K6SD2(KFL6),K7SD2(KFL7),K8SD2(KFL8),K9SD2(KFL9),
     :          KWSD2(6,KFLW),LDELSD2(KFLW,2),SMVRSD2(KFL1)
      DIMENSION K6SE2(KFL6),K7SE2(KFL7),K8SE2(KFL8),K9SE2(KFL9),
     :          KWSE2(6,KFLW),LDELSE2(KFLW,2),SMVRSE2(KFL1)
      DIMENSION K6AD2(KFL6),K7AD2(KFL7),K8AD2(KFL8),K9AD2(KFL9),
     :          KWAD2(6,KFLW),LDELAD2(KFLW,2),SMVRAD2(KFL1)
      DIMENSION K6AE2(KFL6),K7AE2(KFL7),K8AE2(KFL8),K9AE2(KFL9),
     :          KWAE2(6,KFLW),LDELAE2(KFLW,2),SMVRAE2(KFL1)
C
      DIMENSION J6PSD2(KFLV),J7PSD2(KFLV),J8PSD2(KFLV),J9PSD2(KFLV),
     :    JWRDSD2(6,KFLW),NBJSD2(KFLN),NB6JSD2(KFLN),K6CPSD2(KFLN),
     :    K7CPSD2(KFLN),K8CPSD2(KFLN),K9CPSD2(KFLN),JSM6SD2(KFLS),
     :    JSM4SD2(KFLS,KFLW),JSM5SD2(KFLS,KFLW),IN6JSD2(KFLW)
      DIMENSION J6PSE2(KFLV),J7PSE2(KFLV),J8PSE2(KFLV),J9PSE2(KFLV),
     :    JWRDSE2(6,KFLW),NBJSE2(KFLN),NB6JSE2(KFLN),K6CPSE2(KFLN),
     :    K7CPSE2(KFLN),K8CPSE2(KFLN),K9CPSE2(KFLN),JSM6SE2(KFLS),
     :    JSM4SE2(KFLS,KFLW),JSM5SE2(KFLS,KFLW),IN6JSE2(KFLW)
      DIMENSION J6PAD2(KFLV),J7PAD2(KFLV),J8PAD2(KFLV),J9PAD2(KFLV),
     :    JWRDAD2(6,KFLW),NBJAD2(KFLN),NB6JAD2(KFLN),K6CPAD2(KFLN),
     :    K7CPAD2(KFLN),K8CPAD2(KFLN),K9CPAD2(KFLN),JSM6AD2(KFLS),
     :    JSM4AD2(KFLS,KFLW),JSM5AD2(KFLS,KFLW),IN6JAD2(KFLW)
      DIMENSION J6PAE2(KFLV),J7PAE2(KFLV),J8PAE2(KFLV),J9PAE2(KFLV),
     :    JWRDAE2(6,KFLW),NBJAE2(KFLN),NB6JAE2(KFLN),K6CPAE2(KFLN),
     :    K7CPAE2(KFLN),K8CPAE2(KFLN),K9CPAE2(KFLN),JSM6AE2(KFLS),
     :    JSM4AE2(KFLS,KFLW),JSM5AE2(KFLS,KFLW),IN6JAE2(KFLW)
 
      DIMENSION SPINDT(2),SPINEX(2),ANGDIR(2),ANGEX(2)
      DIMENSION BDIRCT(2),BEXCHG(2)
      COMMON/INFORM/IREAD,IWRITE,IOUT,ISC0,ISC1,ISC2,ISC3,JSC(3),IALL
      COMMON/DEBUG/IBUG1,IBUG2,IBUG3,NBUG6,NBUG7,IFULL
      COMMON/KRON/IDEL(10,10)
      COMMON/TERMS/NROWS,ITAB(24),JTAB(24),NTAB(333)
      COMMON/MEDEFN/IHSH,NJ(10),LJ(10),NOSH1(10),NOSH2(10),J1QN1(19,3),
     :     J1QN2(19,3),IJFUL(10)
      COMMON/MVALUE/M1,M2,M3,M4,M5,M6,M7,M8,M9,M10,M11,M12,M13,M14,M15,
     :     M16,M17,M18,M19,M20
      COMMON/XATION/KD1,KD2,KE1,KE2,MULTD,MULTE,
     :     AMULT1(19),AMULT2(19),BMULT1(19),BMULT2(19)
      COMMON/NJLJ/NRHO,LRHO,NSIG,LSIG,NRHOP,LRHOP,NSIGP,LSIGP
      COMMON/INTERM/J1BAR1(10,3),J1BAR2(10,3),J1TLD1(3),J1TLD2(3)
      COMMON/COUPLE/NJ1S,NJ23S,J1(KFL1),J2(KFL2,3),J3(KFL2,3),FREE(KFL1)
 
      COMMON/OVRLAP/MU,NU,MUP,NUP,NONORT,NOVLPS,IROWMU,IROWNU,ICOLMU,
     : ICOLNU,NORTH,IORDER,NCALLS,LMU,LNU,LMUP,LNUP,JMU,JNU,JMUP,JNUP,
     :     IORTH(NORD)
      COMMON/MACHOR/IMATCH(2)
      COMMON/OVRINT/IOVEL(2),IOVER(2),IOVEP(2)
      COMMON/FACT/GAM(100)
*
*
* === DETERMINES POTENTIAL MATRIX ELEMENT FOR GIVEN RHO,SIG,RHOP,SIGP
*
*
  301 FORMAT(21H NO POSSIBLE K-VALUES)
  302 FORMAT(66H SPECTATOR QUANTUM NUMBERS NOT DIAGONAL FOR NON-INTERACT
     :ING SHELLS)
  303 FORMAT(24H DIRECT SPIN INTEGRALS =,2F12.6)
  304 FORMAT(26H EXCHANGE SPIN INTEGRALS =,2F12.6)
  305 FORMAT(23H DIRECT SPIN INTEGRAL =,F10.6)
  306 FORMAT(25H EXCHANGE SPIN INTEGRAL =,F10.6)
  307 FORMAT(6H NJ,LJ,4(I8,I4))
  308 FORMAT(6H KD1 =,I4,6H KD2 =,I4,6H KE1 =,I4,6H KE2 =,I4)
  309 FORMAT(56H ROWS OF TERM TABLE CONTAINING PARENTS ARE, RESPECTIVELY
     :,2(I6,I3))
  310 FORMAT(26H DIRECT ANGULAR INTEGRAL =,F10.6)
  311 FORMAT(27H DIRECT ANGULAR INTEGRALS =,2F12.6)
  312 FORMAT(28H EXCHANGE ANGULAR INTEGRAL =,F10.6)
  313 FORMAT(29H EXCHANGE ANGULAR INTEGRALS =,2F12.6)
*
* --- SET USEFUL CONSTANTS
*
      IORDER=1
      IOTHER=2
      IF(NOVLPS.LT.2) GO TO 960
      IF(MU.EQ.NU.AND.MUP.EQ.NUP) GO TO 960
      NONORT=1
      NKAP = 0
      GO TO 961
  960 NONORT=0
  961 M1=ISIG-IRHO
      M2=ISIGP-IRHOP
      M19=IRHO-IRHOP
      M20=ISIG-ISIGP
      MULTD=0
      MULTE=0
 
      ISPND1=0
      ISPND2=0
      ISPNE1=0
      ISPNE2=0
      IANGD1=0
      IANGD2=0
      IANGE1=0
      IANGE2=0
*
*     JSNDIR,JANGDI=0  IMPLIES APPROPRIATE J2,J3 ARRAYS FOR CALL OF
*     NJSYM HAVE NOT BEEN SET
*
      JSNDIR=0
      JANGDI=0
*
*     Set  the FAIL parameters .FALSE> initially so that NJGRAF can
*     be called
*
      FAILSD1 = .FALSE.
      FAILSE1 = .FALSE.
      FAILAD1 = .FALSE.
      FAILAE1 = .FALSE.
      FAILSD2 = .FALSE.
      FAILSE2 = .FALSE.
      FAILAD2 = .FALSE.
      FAILAE2 = .FALSE.
*
* --- SET N,L VALUES OF INTERACTING SHELLS
*
      NRHO=NJ(IRHO)
      LRHO=LJ(IRHO)
      NSIG=NJ(ISIG)
      LSIG=LJ(ISIG)
      NRHOP=NJ(IRHOP)
      LRHOP=LJ(IRHOP)
      NSIGP=NJ(ISIGP)
      LSIGP=LJ(ISIGP)
      IF(IBUG2-1) 160,161,161
  161 WRITE(IWRITE,307) NRHO,LRHO,NSIG,LSIG,NRHOP,LRHOP,NSIGP,LSIGP
*
* --- EVALUATE MULTIPLICATIVE FACTORS
*
  160 IL=IDEL(IRHO,ISIG)
      IR=IDEL(IRHOP,ISIGP)
      Q=NOSH1(IRHO)*(NOSH1(ISIG)-IL)*NOSH2(IRHOP)*(NOSH2(ISIGP)-IR)
      QD=0.5D0
  901 XMULT=SQRT(Q)*QD
      ADIRCT=DFLOAT(2-IL-IR+IL*IR)/SQRT(DFLOAT(4*LSIG*LRHOP+2*LSIG+2*LR
     :HOP+1))
      IEXCHG=2-IL-IR
      AEXCHG=DFLOAT(IEXCHG)/SQRT(DFLOAT(4*LSIG*LSIGP+2*LSIG+2*LSIGP+1))
  902 DO 13 J=1,IHSH
      NBAR(J)=NOSH1(J)-IDEL(J,IRHO)-IDEL(J,ISIG)
   13 CONTINUE
      IF (MUP .NE. 0) NBAR(MUP) = NOSH2(MUP)
     :   - IDEL(MUP,IRHOP) - IDEL(MUP,ISIGP)
      IF (NUP .NE. 0) NBAR(NUP) = NOSH2(NUP)
     :   - IDEL(NUP,IRHOP) - IDEL(NUP,ISIGP)
      IF(NOVLPS.LT.2) GO TO 930
      IF(MU.EQ.NU) THEN
         NKAP=NBAR(MUP)*NBAR(NUP)
         GO TO 931
      ELSEIF(MUP.EQ.NUP) THEN
         NKAP=NBAR(MU)*NBAR(NU)
         GO TO 931
      ELSE
         NKAP=NBAR(MU)*NBAR(NU)
         GO TO 930
      ENDIF
  931 N1=IOVEP(1)
      N2=IOVEP(2)
      N3=N1+N2
      XMULT=XMULT*EXP((GAM(N3+1)-GAM(N1+1)-GAM(N2+1))/2.D0)
  930 IDELP=0
      IF(M1) 14,15,14
   14 K1=IRHO+1
      DO 16 J=K1,ISIG
      IF (J.EQ.MUP .OR. J.EQ.NUP) GO TO 16
      IDELP=IDELP+NBAR(J)
   16 CONTINUE
   15 IF(M2) 17,18,17
   17 K1=IRHOP+1
      DO 19 J=K1,ISIGP
      IF ( J.EQ.MU .OR. J.EQ.NU) GO TO 19
      IDELP=IDELP+NBAR(J)
   19 CONTINUE
   18 IF(NOVLPS.EQ.0) GO TO 903
      IF(NOVLPS.EQ.1) GO TO 925
      IF (MU.EQ.NU) GO TO 904
      IF (MUP.EQ.NUP) GO TO 905
      IF (NBAR(MU).EQ.NBAR(NUP)) GO TO 906
      MU1N = MIN0(MU,MUP) + 1
      MU1X = MAX0(MU,MUP) - 1
      MU2N = MIN0(NU,NUP) + 1
      MU2X = MAX0(NU,NUP) - 1
      MU1=MU
      NU1=NU
      GO TO 933
  906 MU1N = MIN0(MU,NUP) + 1
      MU1X = MAX0(MU,NUP) - 1
      MU2N = MIN0(NU,MUP) + 1
      MU2X = MAX0(NU,MUP) - 1
      MU1=MU
      NU1=NU
      GO TO 933
  905 MU1N = MIN0(MU,MUP) + 1
      MU1X = MAX0(MU,MUP) - 1
      MU2N = MIN0(MUP,NU) + 1
      MU2X = MAX0(MUP,NU) - 1
      MU1=MU
      NU1=NU
      GO TO 933
  904 MU1N = MIN0(MU,MUP) + 1
      MU1X = MAX0(MU,MUP) - 1
      MU2N = MIN0(MU,NUP) + 1
      MU2X = MAX0(MU,NUP) - 1
      MU1=MUP
      NU1=NUP
  933 IF (MU1N .GT. MU1X) GO TO 934
      DO 962 J = MU1N,MU1X
      IF (J.EQ.MU .OR. J.EQ.NU .OR. J.EQ.MUP .OR. J.EQ.NUP) GO TO 962
      IDELP = IDELP + NBAR(J)*NBAR(MU1)
  962 CONTINUE
  934 IF (MU2N .GT. MU2X) GO TO 903
      DO 963 J = MU2N,MU2X
      IF (J.EQ.MU .OR. J.EQ.NU .OR. J.EQ.MUP .OR. J.EQ.NUP) GO TO 963
      IDELP = IDELP + NBAR(J)*NBAR(NU1)
  963 CONTINUE
      GO TO 903
  925 MUMIN1=MIN0(MU,MUP)+1
      MUMAX1=MAX0(MU,MUP)-1
      IF(MUMIN1.GT.MUMAX1) GO TO 903
      DO 927 J=MUMIN1,MUMAX1
      IDELP=IDELP+NBAR(J)*NBAR(MU)
  927 CONTINUE
  903 XMULT=XMULT*DFLOAT( (-1)**IDELP )
*
* --- DETERMINES RANGES OF K FOR DIRECT AND EXCHANGE INTEGRALS
*     TRIANGULAR RELATIONS LIMIT (K+1) VALUES TO LIE BETWEEN KD1 AND KD2
*     FOR =DIRECT= INTEGRALS,  KE1 AND KE2 FOR =EXCHANGE= INTEGRALS
*
      K1=IABS(LSIG-LSIGP)
      K2=LSIG+LSIGP
      K3=IABS(LRHO-LRHOP)
      K4=LRHO+LRHOP
      KD1=MAX0(K1,K3)+1
      KD2=MIN0(K2,K4)+1
      K1=IABS(LRHOP-LSIG)
      K2=LRHOP+LSIG
      K3=IABS(LRHO-LSIGP)
      K4=LRHO+LSIGP
      KE1=MAX0(K1,K3)+1
      KE2=MIN0(K2,K4)+1
      IF(IBUG2-1) 612,613,613
  613 WRITE(IWRITE,308) KD1,KD2,KE1,KE2
  612 IF(KD1-KD2) 213,213,211
  211 IF(KE1-KE2) 213,213,212
  212 IF(IBUG2-1) 400,401,401
  401 WRITE(IWRITE,301)
  400 RETURN
*
* --- ZEROIZE MULTIPLYING FACTORS FOR ALLOWED K-VALUES.  THE LOWEST
*     VALUES  KD1  AND  KD2  ARE ALWAYS ALLOWED (UNLESS THEY ARE
*     GREATER THAN KD2 AND KE2 RESPECTIVELY).  ALLOWED K-VALUES THEN
*     STEP BY 2 TO SATISFY THE EVEN CONDITION OF THE REDUCED MATRIX
*     ELEMENTS, WHICH ARE THEN CALCULATED AND STORED
*
  213 IF(KD1-KD2) 231,231,232
  231 DO 230 JK1=KD1,KD2,2
      K=JK1-1
      AMULT1(JK1)=0.D0
      AMULT2(JK1)=0.D0
      RMEDIR(JK1)=RME(LRHO,LRHOP,K)*RME(LSIG,LSIGP,K)
  230 CONTINUE
  232 IF(KE1-KE2) 233,233,241
  233 DO 234 JK1=KE1,KE2,2
      K=JK1-1
      BMULT1(JK1)=0.D0
      BMULT2(JK1)=0.D0
      RMEEX(JK1)=RME(LRHO,LSIGP,K)*RME(LSIG,LRHOP,K)
  234 CONTINUE
*
* --- SET SENIORITY, S AND L VALUES OF SPECTATOR SHELLS
*
  241 DO 26 J=1,IHSH
      IF(IRHO-J) 27,29,27
   27 IF(ISIG-J) 28,29,28
   28 DO 128 K=1,3
      J1BAR1(J,K)=J1QN1(J,K)
  128 CONTINUE
   29 IF(IRHOP-J) 30,26,30
   30 IF(ISIGP-J) 31,26,31
   31 DO 181 K=1,3
      J1BAR2(J,K)=J1QN2(J,K)
  181 CONTINUE
      IF (J.EQ.IRHO .OR. J.EQ.ISIG .OR. J.EQ.MU .OR. J.EQ.NU
     :    .OR. J.EQ.MUP .OR. J.EQ.NUP ) GO TO 26
*
*     CHECK THAT SPECTATOR SHELLS HAVE SAME RESPECTIVE QUANTUM NUMBERS
*
   33 DO 183 K=1,3
      IF(J1BAR1(J,K)-J1BAR2(J,K)) 402,183,402
  183 CONTINUE
   26 CONTINUE
      GO TO 405
  402 IF(IBUG2-1) 404,403,403
  403 WRITE(IWRITE,302)
  404 RETURN
*
* --- FROM WHICH ROWS OF NTAB DO WE FIND THE QUANTUM NUMBERS WITH BARS
*     OR TILDES
*
  405 NELCTS=NOSH1(ISIG)
      K2=NTAB1(NELCTS,LSIG+1)
  361 IF(M1) 20,21,20
   21 NELCTS=NOSH1(ISIG)-1
      GO TO 22
   20 NELCTS=NOSH1(IRHO)
   22 K1=NTAB1(NELCTS,LRHO+1)
  363 NELCTS=NOSH2(ISIGP)
      K4=NTAB1(NELCTS,LSIGP+1)
  365 IF(M2) 23,24,23
   24 NELCTS=NOSH2(ISIGP)-1
      GO TO 25
   23 NELCTS=NOSH2(IRHOP)
   25 K3=NTAB1(NELCTS,LRHOP+1)
  367 IF(IBUG2-1) 59,49,49
   49 WRITE(IWRITE,309) K1,K2,K3,K4
   59 KK1=ITAB(K1)
      KK2=ITAB(K2)
      KK3=ITAB(K3)
      KK4=ITAB(K4)
*
* === SUM OVER QUANTUM NUMBERS WITH BARS OR TILDES
*
      DO 151 JJ2=1,KK2
*
* --- TEST TO SEE WHICH PARENT TERMS ARE ALLOWABLE.  ONLY TEST THIS ON
*     L  AND  S  VALUES AT THIS STAGE, BY MEANS OF TRIANGULAR CONDITIONS
*     FOR TWICE THE QUANTUM NUMBERS, IN ORDER TO USE ONLY INTEGER
*     QUANTITIES
*
      IN3=2*LSIG
      IJK2=3*(JJ2-1)+JTAB(K2)
      DO 131 K=2,3
      IN1=NTAB(IJK2+K)-1
      IN2=J1QN1(ISIG,K)-1
      IF(IN1-IN2-IN3) 130,130,151
  130 IF(IN1-IABS(IN2-IN3)) 151,140,140
  140 IN3=1
  131 CONTINUE
  351 DO 152 JJ1=1,KK1
      IN3=2*LRHO
      IJK1=3*(JJ1-1)+JTAB(K1)
  162 DO 132 K=2,3
      IN1=NTAB(IJK1+K)-1
      IF(M1) 141,142,141
  141 IN2=J1QN1(IRHO,K)-1
      GO TO 143
  142 IN2=NTAB(IJK2+K)-1
  143 IF(IN1-IN2-IN3) 144,144,152
  144 IF(IN1-IABS(IN2-IN3)) 152,145,145
  145 IN3=1
  132 CONTINUE
  352 DO 153 JJ4=1,KK4
      IN3=2*LSIGP
      IJK4=3*(JJ4-1)+JTAB(K4)
      DO 133 K=2,3
      IN1=NTAB(IJK4+K)-1
      IN2=J1QN2(ISIGP,K)-1
      IF(IN1-IN2-IN3) 146,146,153
  146 IF(IN1-IABS(IN2-IN3)) 153,147,147
  147 IN3=1
  133 CONTINUE
  353 DO 154 JJ3=1,KK3
      I1BRNO = 1
      I2BRNO = 1
      IN3=2*LRHOP
      IJK3=3*(JJ3-1)+JTAB(K3)
  137 DO 134 K=2,3
      IN1=NTAB(IJK3+K)-1
      IF(M2) 138,139,138
  138 IN2=J1QN2(IRHOP,K)-1
      GO TO 148
  139 IN2=NTAB(IJK4+K)-1
  148 IF(IN1-IN2-IN3) 149,149,154
  149 IF(IN1-IABS(IN2-IN3)) 154,150,150
  150 IN3=1
  134 CONTINUE
*
*     SUMMATIONS NOW PERFORMED OVER ALLOWED QUANTUM NUMBERS
*     THE TILDES CORRESPOND TO  IRHO=ISIG  AND/OR IRHOP=ISIGP
*
* --- SET THE REMAINING QUANTUM NUMBERS WITH BARS OR TILDES
*
  354 DO 35 K=1,3
      J1BAR1(IRHO,K)=NTAB(IJK1+K)
      J1BAR2(IRHOP,K)=NTAB(IJK3+K)
      IF(M1) 36,37,36
   36 J1BAR1(ISIG,K)=NTAB(IJK2+K)
      GO TO 38
   37 J1TLD1(K)=NTAB(IJK2+K)
   38 IF( M2) 39,40,39
   39 J1BAR2(ISIGP,K)=NTAB(IJK4+K)
      GO TO 35
   40 J1TLD2(K)=NTAB(IJK4+K)
   35 CONTINUE
*
* --- IS POTENTIAL DIAG. IN BARRED QU. NOS. FOR INTERACTING SHELLS
*
      I5=0
      I=ISIG
      GO TO 50
   42 I=IRHO
      IF( M1) 43,44,43
   43 GO TO 50
   44 I5=I5+1
   45 I=ISIGP
      GO TO 50
   46 I=IRHOP
      IF(M2) 47,48,47
   47 GO TO 50
   50 I5=I5+1
      IF(I.EQ.MU.OR.I.EQ.NU.OR.I.EQ.MUP.OR.I.EQ.NUP) GO TO 360
      DO 51 K=1,3
      IF(J1BAR1(I,K)-J1BAR2(I,K)) 154,51,154
   51 CONTINUE
  360 GO TO (42,45,46,48),I5
   48 PICFP=1.D0
*
* --- EVALUATE FRACTIONAL PARENTAGE COEFFICIENTS
*
      IF (NOVLPS .NE. 1) GO TO 935
      DO 936 K = 1,3
      IF (J1BAR1(MU,K) .NE. J1BAR2(MUP,K)) GO TO 154
  936 CONTINUE
935   I=1
      CALL MUMDAD (I,ISIG,IRHO,M1,CFPLHS)
      PICFP=PICFP*CFPLHS
      IF(DABS(PICFP).LT.1.D-14) GO TO 154
   53 I=2
      CALL MUMDAD(I,ISIGP,IRHOP,M2,CFPRHS)
      PICFP=PICFP*CFPRHS
      IF(DABS(PICFP).LT.1.D-14) GO TO 154
*
* === SET UP J1,J2,J3 AND EVALUATE RECOUPLING COEFFICIENTS
*
* --- FIRST OF ALL, THE SPIN COEFFICIENTS
*
   55 I=3
      CALL SETJ1(I,IRHO,ISIG,IRHOP,ISIGP,ISPND1,ISPND2,ISPNE1,ISPNE2,
     :           KK1,KK2,KK3,KK4)
  269 CALL J23SPN(IRHO,ISIG,IRHOP,ISIGP,JSNDIR)
*
* --- DIRECT SPIN INTEGRAL
*
  570 IF(KD1-KD2) 89,89,90
   90 SPINDT(1)=0.D0
      IF(NOVLPS.LT.2) SPINDT(2)=0.D0
      GO TO 78
   89 IF (NOVLPS .EQ. 1) GO TO 978
      IF(NOVLPS.LT.2 .OR. MU.EQ.NU .OR. MUP.EQ.NUP) GO TO 950
      IF(LMU.NE.LMUP.OR.LNU.NE.LNUP) GO TO 947
      DO 937 K = 1,3
      IF (J1BAR1(MU,K) .NE. J1BAR2(MUP,K)) GO TO 939
      IF (J1BAR1(NU,K) .NE. J1BAR2(NUP,K)) GO TO 939
  937 CONTINUE
      GO TO 950
  939 I1BRNO = 0
      GO TO 947
  978 DO 979 K = 1,3
      IF (J1BAR1(MU,K) .NE. J1BAR2(MUP,K)) GO TO 939
  979 CONTINUE
 
  950 CONTINUE
 
      IF (.NOT.FAILSD1) THEN
        IF (ISPND1.EQ.0) THEN
CDBG      print*,' in failsd1  nj1s = ',nj1s,' nj23 = ',nj23s
CDBG      print*,' j1 = ',j1
CDBG      print*,' j2 = ',j2
CDBG      print*,' j3 = ',j3
          CALL NJGRAF(SPINDT(1),FAILSD1)
          ISPND1=1
          IF(FAILSD1) GO TO 947
          CALL KNJ(J6CSD1,J7CSD1,J8CSD1,J9CSD1,JWCSD1,K6SD1,K7SD1,
     :             K8SD1,K9SD1,KWSD1,JDELSD1,LDELSD1,SMVRSD1,MPSD1,
     :             J6PSD1,J7PSD1,J8PSD1,J9PSD1,JWRDSD1,NLSMSD1,NBJSD1,
     :             NB6JSD1,K6CPSD1,K7CPSD1,K8CPSD1,K9CPSD1,JSM4SD1,
     :             JSM5SD1,JSM6SD1,IN6JSD1)
        ELSE
          CALL GENSUM(J6CSD1,J7CSD1,J8CSD1,J9CSD1,JWCSD1,K6SD1,K7SD1,
     :             K8SD1,K9SD1,KWSD1,JDELSD1,LDELSD1,SMVRSD1,MPSD1,
     :             J6PSD1,J7PSD1,J8PSD1,J9PSD1,JWRDSD1,NLSMSD1,NBJSD1,
     :             NB6JSD1,K6CPSD1,K7CPSD1,K8CPSD1,K9CPSD1,JSM4SD1,
     :             JSM5SD1,JSM6SD1,IN6JSD1,SPINDT(1))
C
        ENDIF
        GO TO 948
      ENDIF
  947 SPINDT(1) = 0.D0
  948 IF(NOVLPS.LT.2) GO TO 78
      IF(MU.EQ.NU.OR.MUP.EQ.NUP) GO TO 940
      IF(LMU.NE.LNUP.OR.LNU.NE.LMUP) GO TO 949
      DO 938 K = 1,3
      IF (J1BAR1(MU,K) .NE. J1BAR2(NUP,K)) GO TO 940
      IF (J1BAR1(NU,K) .NE. J1BAR2(MUP,K)) GO TO 940
  938 CONTINUE
      GO TO 941
  940 I2BRNO = 0
      GO TO 949
  941 I=3
      CALL MODJ23(I)
 
      IF (.NOT.FAILSD2) THEN
        IF (ISPND2.EQ.0) THEN
CDBG      print*,' in failsd2  nj1s = ',nj1s,' nj23 = ',nj23s
CDBG      print*,' j1 = ',j1
CDBG      print*,' j2 = ',j2
CDBG      print*,' j3 = ',j3
          CALL NJGRAF(SPINDT(2),FAILSD2)
          ISPND2=1
          IF(FAILSD2) GO TO 949
          CALL KNJ(J6CSD2,J7CSD2,J8CSD2,J9CSD2,JWCSD2,K6SD2,K7SD2,
     :             K8SD2,K9SD2,KWSD2,JDELSD2,LDELSD2,SMVRSD2,MPSD2,
     :             J6PSD2,J7PSD2,J8PSD2,J9PSD2,JWRDSD2,NLSMSD2,NBJSD2,
     :             NB6JSD2,K6CPSD2,K7CPSD2,K8CPSD2,K9CPSD2,JSM4SD2,
     :             JSM5SD2,JSM6SD2,IN6JSD2)
        ELSE
          CALL GENSUM(J6CSD2,J7CSD2,J8CSD2,J9CSD2,JWCSD2,K6SD2,K7SD2,
     :             K8SD2,K9SD2,KWSD2,JDELSD2,LDELSD2,SMVRSD2,MPSD2,
     :             J6PSD2,J7PSD2,J8PSD2,J9PSD2,JWRDSD2,NLSMSD2,NBJSD2,
     :             NB6JSD2,K6CPSD2,K7CPSD2,K8CPSD2,K9CPSD2,JSM4SD2,
     :             JSM5SD2,JSM6SD2,IN6JSD2,SPINDT(2))
C
        ENDIF
        GO TO 78
      ENDIF
  949 SPINDT(2) = 0.D0
   78 IF(IBUG2-1) 91,170,170
  170 IF(NOVLPS.LT.2) GO TO 907
      WRITE(IWRITE,303) SPINDT(1),SPINDT(2)
      GO TO 91
  907 WRITE(IWRITE,305) SPINDT(1)
*
*     IEXCHG IS ZERO WHENEVER  M1=0=M2 , IN WHICH CASE THE EXCHANGE
*     INTEGRAL HAS ZERO COEFFICIENT.  THERE IS THEN NO POINT IN
*     CALCULATING THIS INTEGRAL, AND SPINEX IS SET ZERO (AT STATEMENT
*     93) AS A MARKER OF THIS SITUATION
*
   91 IF(IEXCHG.EQ.0) GO TO 93
*
* --- MODIFY J2 AND J3 TO CALCULATE THE EXCHANGE SPIN INTEGRAL
*
      I=1
      CALL MODJ23(I)
*
* --- EXCHANGE SPIN INTEGRAL
*
      IF(KE1-KE2) 92,92,93
   93 SPINEX(1)=0.D0
      IF(NOVLPS.LT.2) SPINEX(2)=0.D0
      GO TO 94
   92 IF (I1BRNO .EQ. 0) GO TO 953
      IF(NOVLPS.LT.2) GO TO 954
      IF(MU.EQ.NU.OR.MUP.EQ.NUP) GO TO 954
      IF(IMATCH(1).EQ.0) GO TO 953
  951 IF(LMU.NE.LMUP.OR.LNU.NE.LNUP.OR.I1BRNO.EQ.0) GO TO 953
 
  954 CONTINUE
      IF (.NOT.FAILSE1) THEN
        IF (ISPNE1.EQ.0) THEN
CDBG      print*,' in failse1  nj1s = ',nj1s,' nj23 = ',nj23s
CDBG      print*,' j1 = ',j1
CDBG      print*,' j2 = ',j2
CDBG      print*,' j3 = ',j3
          CALL NJGRAF(SPINEX(1),FAILSE1)
          ISPNE1=1
          IF(FAILSE1) GO TO 953
          CALL KNJ(J6CSE1,J7CSE1,J8CSE1,J9CSE1,JWCSE1,K6SE1,K7SE1,
     :             K8SE1,K9SE1,KWSE1,JDELSE1,LDELSE1,SMVRSE1,MPSE1,
     :             J6PSE1,J7PSE1,J8PSE1,J9PSE1,JWRDSE1,NLSMSE1,NBJSE1,
     :             NB6JSE1,K6CPSE1,K7CPSE1,K8CPSE1,K9CPSE1,JSM4SE1,
     :             JSM5SE1,JSM6SE1,IN6JSE1)
        ELSE
          CALL GENSUM(J6CSE1,J7CSE1,J8CSE1,J9CSE1,JWCSE1,K6SE1,K7SE1,
     :             K8SE1,K9SE1,KWSE1,JDELSE1,LDELSE1,SMVRSE1,MPSE1,
     :             J6PSE1,J7PSE1,J8PSE1,J9PSE1,JWRDSE1,NLSMSE1,NBJSE1,
     :             NB6JSE1,K6CPSE1,K7CPSE1,K8CPSE1,K9CPSE1,JSM4SE1,
     :             JSM5SE1,JSM6SE1,IN6JSE1,SPINEX(1))
C
        ENDIF
        GO TO 955
      ENDIF
  953 SPINEX(1)=0.D0
  955 IF(NOVLPS.LT.2) GO TO 94
      IF(IMATCH(2).EQ.0) GO TO 958
  956 IF(LMU.NE.LNUP.OR.LNU.NE.LMUP.OR.I2BRNO.EQ.0) GO TO 958
  959 I=5
      CALL MODJ23(I)
 
      IF (.NOT.FAILSE2) THEN
        IF (ISPNE2.EQ.0) THEN
CDBG      print*,' in failse2  nj1s = ',nj1s,' nj23 = ',nj23s
CDBG      print*,' j1 = ',j1
CDBG      print*,' j2 = ',j2
CDBG      print*,' j3 = ',j3
          CALL NJGRAF(SPINEX(2),FAILSE2)
          ISPNE2=1
          IF(FAILSE2) GO TO 958
          CALL KNJ(J6CSE2,J7CSE2,J8CSE2,J9CSE2,JWCSE2,K6SE2,K7SE2,
     :             K8SE2,K9SE2,KWSE2,JDELSE2,LDELSE2,SMVRSE2,MPSE2,
     :             J6PSE2,J7PSE2,J8PSE2,J9PSE2,JWRDSE2,NLSMSE2,NBJSE2,
     :             NB6JSE2,K6CPSE2,K7CPSE2,K8CPSE2,K9CPSE2,JSM4SE2,
     :             JSM5SE2,JSM6SE2,IN6JSE2)
        ELSE
          CALL GENSUM(J6CSE2,J7CSE2,J8CSE2,J9CSE2,JWCSE2,K6SE2,K7SE2,
     :             K8SE2,K9SE2,KWSE2,JDELSE2,LDELSE2,SMVRSE2,MPSE2,
     :             J6PSE2,J7PSE2,J8PSE2,J9PSE2,JWRDSE2,NLSMSE2,NBJSE2,
     :             NB6JSE2,K6CPSE2,K7CPSE2,K8CPSE2,K9CPSE2,JSM4SE2,
     :             JSM5SE2,JSM6SE2,IN6JSE2,SPINEX(2))
C
        ENDIF
        GO TO 94
      ENDIF
  958 SPINEX(2) = 0.D0
   94 IF(IBUG2-1) 171,172,172
  172 IF(NOVLPS.LT.2) GO TO 908
      WRITE(IWRITE,304) SPINEX(1),SPINEX(2)
      GO TO 171
  908 WRITE(IWRITE,306) SPINEX(1)
*
* --- MULTIPLY SPIN INTEGRALS BY PRODUCT OF FRACTIONAL PARENTAGE COEFFS
*
  171 BDIRCT(1)=SPINDT(1)*PICFP
      BEXCHG(1)=SPINEX(1)*PICFP
      IF(NOVLPS.LT.2) GO TO 916
      BDIRCT(2)=SPINDT(2)*PICFP
      BEXCHG(2)=SPINEX(2)*PICFP
*
* --- THE ANGULAR RECOUPLING COEFFICIENTS
*     SET J1,J2,J3   (COMPARE SPIN INTEGRALS)
*
*     IF BOTH SPIN INTEGRALS ARE ZERO,  THERE IS NO PURPOSE IN
*     CALCULATING THE ANGULAR INTEGRALS
*
  909 DT=DABS(SPINDT(1))+DABS(SPINDT(2))
      EX=DABS(SPINEX(1))+DABS(SPINEX(2))
      GO TO 917
  916 DT=DABS(SPINDT(1))
      EX=DABS(SPINEX(1))
  917 DPE=DT+EX
      IF(DPE.LT.1.0D-14) GO TO 154
*
   87 I=2
      CALL SETJ1(I,IRHO,ISIG,IRHOP,ISIGP,IANGD1,IANGD2,IANGE1,IANGE2,
     :           KK1,KK2,KK3,KK4)
      IF (IANGD1.EQ.0) CALL J23ANG(IRHO,ISIG,IRHOP,ISIGP,JANGDI)
*
*     IF THE DIRECT SPIN RECOUPLING COEFFICIENT IS ZERO, WE NEED NOT
*     CALCULATE THE CORRESPONDING ORBITAL RECOUPLING COEFFICIENT
*
      IF(DT.LT.1.0D-14) GO TO 121
*
* --- DIRECT ANGULAR INTEGRAL
*
*     CONSIDER ALL ALLOWED K-VALUES
*
      IF(KD2.NE.KD1) THEN
          FREE(NJ1S)=.TRUE.
      ELSE
          FREE(NJ1S)=.FALSE.
      ENDIF
*
      DO 114 JK1=KD1,KD2,2
          J1(NJ1S)=2*JK1-1
          IF(.NOT.FAILAD1.AND.DABS(SPINDT(1)).GT.1.0D-14) THEN
              IF(IANGD1.EQ.0) THEN
                  CALL J23ANG(IRHO,ISIG,IRHOP,ISIGP,JANGDI)
CDBG  print*,' in failad1 (do 114)  nj1s = ',nj1s,' nj23 = ',nj23s
CDBG      print*,' j1 = ',j1
CDBG      print*,' j2 = ',j2
CDBG      print*,' j3 = ',j3
                  CALL NJGRAF(ANGDIR(1),FAILAD1)
                  IANGD1=1
                  IF(.NOT.FAILAD1) THEN
                  CALL KNJ(J6CAD1,J7CAD1,J8CAD1,J9CAD1,JWCAD1,
     :                     K6AD1,K7AD1,K8AD1,K9AD1,KWAD1,JDELAD1,
     :                     LDELAD1,SMVRAD1,MPAD1,J6PAD1,J7PAD1,
     :                     J8PAD1,J9PAD1,JWRDAD1,NLSMAD1,NBJAD1,
     :                     NB6JAD1,K6CPAD1,K7CPAD1,K8CPAD1,K9CPAD1,
     :                     JSM4AD1,JSM5AD1,JSM6AD1,IN6JAD1)
                  ENDIF
              ELSE
                  CALL GENSUM(J6CAD1,J7CAD1,J8CAD1,J9CAD1,JWCAD1,
     :                     K6AD1,K7AD1,K8AD1,K9AD1,KWAD1,JDELAD1,
     :                     LDELAD1,SMVRAD1,MPAD1,J6PAD1,J7PAD1,
     :                     J8PAD1,J9PAD1,JWRDAD1,NLSMAD1,NBJAD1,
     :                     NB6JAD1,K6CPAD1,K7CPAD1,K8CPAD1,K9CPAD1,
     :                     JSM4AD1,JSM5AD1,JSM6AD1,IN6JAD1,ANGDIR(1))
              ENDIF
              AMULT1(JK1)=AMULT1(JK1)+ANGDIR(1)*BDIRCT(1)
              MULTD=1
          ELSE
              ANGDIR(1)=0.0D0
          ENDIF
C
          IF(NOVLPS.LT.2) GO TO 912
C
          IF(.NOT.FAILAD2.AND.DABS(SPINDT(2)).GT.1.0D-14) THEN
              IF(IANGD2.EQ.0) THEN
                  I=4
                  CALL MODJ23(I)
CDBG  print*,' in failad2 (in do 114) nj1s = ',nj1s,' nj23 = ',nj23s
CDBG      print*,' j1 = ',j1
CDBG      print*,' j2 = ',j2
CDBG      print*,' j3 = ',j3
                  CALL NJGRAF(ANGDIR(2),FAILAD2)
                  IANGD2=1
                  IF(.NOT.FAILAD2) THEN
                  CALL KNJ(J6CAD2,J7CAD2,J8CAD2,J9CAD2,JWCAD2,
     :                     K6AD2,K7AD2,K8AD2,K9AD2,KWAD2,JDELAD2,
     :                     LDELAD2,SMVRAD2,MPAD2,J6PAD2,J7PAD2,
     :                     J8PAD2,J9PAD2,JWRDAD2,NLSMAD2,NBJAD2,
     :                     NB6JAD2,K6CPAD2,K7CPAD2,K8CPAD2,K9CPAD2,
     :                     JSM4AD2,JSM5AD2,JSM6AD2,IN6JAD2)
                  ENDIF
              ELSE
                  CALL GENSUM(J6CAD2,J7CAD2,J8CAD2,J9CAD2,JWCAD2,
     :                     K6AD2,K7AD2,K8AD2,K9AD2,KWAD2,JDELAD2,
     :                     LDELAD2,SMVRAD2,MPAD2,J6PAD2,J7PAD2,
     :                     J8PAD2,J9PAD2,JWRDAD2,NLSMAD2,NBJAD2,
     :                     NB6JAD2,K6CPAD2,K7CPAD2,K8CPAD2,K9CPAD2,
     :                     JSM4AD2,JSM5AD2,JSM6AD2,IN6JAD2,ANGDIR(2))
              ENDIF
              AMULT2(JK1)=AMULT2(JK1)+ANGDIR(2)*BDIRCT(2)*(-1)**NKAP
              MULTD=1
          ELSE
              ANGDIR(2)=0.0D0
          ENDIF
  912     IF(IBUG2.LT.1) GO TO 114
          IF(NOVLPS.LT.2) THEN
              WRITE(IWRITE,310) ANGDIR(1)
          ELSE
              WRITE(IWRITE,311) ANGDIR(1),ANGDIR(2)
          ENDIF
  114 CONTINUE
 
*
*     IF THE EXCHANGE SPIN RECOUPLING COEFFICIENT IS ZERO, WE NEED NOT
*     CALCULATE THE CORRESPONDING ORBITAL RECOUPLING COEFFICIENT
*
  121 IF(EX.LT.1.D-14) GO TO 154
*
* --- EXCHANGE ANGULAR INTEGRAL
*
*     CONSIDER ALL ALLOWED K-VALUES
*
*
      IF(KE2.NE.KE1) THEN
          FREE(NJ1S)=.TRUE.
      ELSE
          FREE(NJ1S)=.FALSE.
      ENDIF
C
      DO 115 JK1=KE1,KE2,2
          J1(NJ1S)=2*JK1-1
          IF(.NOT.FAILAE1.AND.DABS(SPINEX(1)).GT.1.0D-14) THEN
              IF(IANGE1.EQ.0) THEN
                  I = 2
                  CALL MODJ23(I)
CDBG  print*,' in failae1 (do 115) nj1s = ',nj1s,' nj23 = ',nj23s
CDBG      print*,' j1 = ',j1
CDBG      print*,' j2 = ',j2
CDBG      print*,' j3 = ',j3
                  CALL NJGRAF(ANGEX(1),FAILAE1)
                  IANGE1=1
                  IF(.NOT.FAILAE1) THEN
                  CALL KNJ(J6CAE1,J7CAE1,J8CAE1,J9CAE1,JWCAE1,
     :                     K6AE1,K7AE1,K8AE1,K9AE1,KWAE1,JDELAE1,
     :                     LDELAE1,SMVRAE1,MPAE1,J6PAE1,J7PAE1,
     :                     J8PAE1,J9PAE1,JWRDAE1,NLSMAE1,NBJAE1,
     :                     NB6JAE1,K6CPAE1,K7CPAE1,K8CPAE1,K9CPAE1,
     :                     JSM4AE1,JSM5AE1,JSM6AE1,IN6JAE1)
                  ENDIF
              ELSE
                  CALL GENSUM(J6CAE1,J7CAE1,J8CAE1,J9CAE1,JWCAE1,
     :                     K6AE1,K7AE1,K8AE1,K9AE1,KWAE1,JDELAE1,
     :                     LDELAE1,SMVRAE1,MPAE1,J6PAE1,J7PAE1,
     :                     J8PAE1,J9PAE1,JWRDAE1,NLSMAE1,NBJAE1,
     :                     NB6JAE1,K6CPAE1,K7CPAE1,K8CPAE1,K9CPAE1,
     :                     JSM4AE1,JSM5AE1,JSM6AE1,IN6JAE1,ANGEX(1))
              ENDIF
              BMULT1(JK1)=BMULT1(JK1)-ANGEX(1)*BEXCHG(1)
              MULTE=1
          ELSE
              ANGEX(1)=0.0D0
          ENDIF
C
          IF(NOVLPS.LT.2) GO TO 919
C
          IF(.NOT.FAILAE2.AND.DABS(SPINEX(2)).GT.1.0D-14) THEN
              IF(IANGE2.EQ.0) THEN
                  I=6
                  CALL MODJ23(I)
CDBG  print*,' in failae2 (do 115)  nj1s = ',nj1s,' nj23 = ',nj23s
CDBG      print*,' j1 = ',j1
CDBG      print*,' j2 = ',j2
CDBG      print*,' j3 = ',j3
                  CALL NJGRAF(ANGEX(2),FAILAE2)
                  IANGE2=1
                  IF(.NOT.FAILAE2) THEN
                  CALL KNJ(J6CAE2,J7CAE2,J8CAE2,J9CAE2,JWCAE2,
     :                     K6AE2,K7AE2,K8AE2,K9AE2,KWAE2,JDELAE2,
     :                     LDELAE2,SMVRAE2,MPAE2,J6PAE2,J7PAE2,
     :                     J8PAE2,J9PAE2,JWRDAE2,NLSMAE2,NBJAE2,
     :                     NB6JAE2,K6CPAE2,K7CPAE2,K8CPAE2,K9CPAE2,
     :                     JSM4AE2,JSM5AE2,JSM6AE2,IN6JAE2)
                  ENDIF
              ELSE
                  CALL GENSUM(J6CAE2,J7CAE2,J8CAE2,J9CAE2,JWCAE2,
     :                     K6AE2,K7AE2,K8AE2,K9AE2,KWAE2,JDELAE2,
     :                     LDELAE2,SMVRAE2,MPAE2,J6PAE2,J7PAE2,
     :                     J8PAE2,J9PAE2,JWRDAE2,NLSMAE2,NBJAE2,
     :                     NB6JAE2,K6CPAE2,K7CPAE2,K8CPAE2,K9CPAE2,
     :                     JSM4AE2,JSM5AE2,JSM6AE2,IN6JAE2,ANGEX(2))
              ENDIF
              BMULT2(JK1)=BMULT2(JK1)-ANGEX(2)*BEXCHG(2)*(-1)**NKAP
              MULTE=1
          ELSE
              ANGEX(2)=0.0D0
          ENDIF
  919     IF(IBUG2.LT.1) GO TO 115
          IF(NOVLPS.LT.2) THEN
              WRITE(IWRITE,312) ANGEX(1)
          ELSE
              WRITE(IWRITE,313) ANGEX(1),ANGEX(2)
          ENDIF
  115 CONTINUE
 
  154 CONTINUE
  153 CONTINUE
  152 CONTINUE
  151 CONTINUE
*
* === INCLUDE MULTIPLICATIVE FACTORS COMMON TO ALL TERMS WITHIN THIS
*     FOUR-FOLD SUMMATION
*
      IF(MULTD) 524,525,524
  524 X1=XMULT*ADIRCT
      DO 518 JK1=KD1,KD2,2
      X2=X1*RMEDIR(JK1)
      AMULT1(JK1)=AMULT1(JK1)*X2
      IF(NONORT.EQ.0) GO TO 518
      AMULT2(JK1)=AMULT2(JK1)*X2
  518 CONTINUE
  525 IF(MULTE) 526,527,526
  526 X1=XMULT*AEXCHG
      DO 519 JK1=KE1,KE2,2
      X2=X1*RMEEX(JK1)
      BMULT1(JK1)=BMULT1(JK1)*X2
      IF(NOVLPS.LT.2) GO TO 519
      BMULT2(JK1)=BMULT2(JK1)*X2
  519 CONTINUE
*
* --- PRINT OUT THE VALUES OF THE COEFFICIENTS OF THE SLATER INTEGRALS
*
*     THE SUBROUTINE PRNTWT IS CALLED FROM RKWTS
*
  527 RETURN
*
* *** DEFINITION OF DIMENSION LIST
*
*     RMEDIR(K),K=KD1,KD2,2  -  DIRECT REDUCED MATRIX ELEMENT PRODUCT
*     RMEEX(K),K=KE1,KE2,2  -  EXCHANGE REDUCED MATRIX ELEMENT PRODUCT
*     KD1,KE1   ARE ALWAYS .GE. 1
*     KD2,KE2   ARE .LE.  1+2*MAX(L-VALUE)
*     NBAR(I), I=1,IHSH  -  NUMBER OF SPECTATOR ELECTRONS IN EACH SHELL
*      THE K6,K7,K8,KW ARRAYS ARE DEFINED IN NJSYM
*
*
      END
*    MCHF_NONH  (Part 1b of 3 -- 1a, 1b and 2) 
*
*     ------------------------------------------------------------------
*       H 0 W T S
*     ------------------------------------------------------------------
*
      SUBROUTINE H0WTS
      IMPLICIT DOUBLE PRECISION (A-H,O-Z)
      PARAMETER (NWD=30,NCD=100,NORD=NWD*(NWD-1)/2)
*
      DIMENSION NBAR(10),IPAR(3),JBAR(10,3),JPBAR(10,3)
      COMMON/DIAGNL/IDIAG,JA,JB
      COMMON /CLOSED/B1ELC(4),NCLOSD,IAJCLD(NWD),LJCLSD(NWD),IBK
      PARAMETER(KFL1=60,KFL2=12)
      LOGICAL FREE, FAIL
      COMMON/COUPLE/NJ1S,NJ23S,J1(KFL1),J2(KFL2,3),J3(KFL2,3),FREE(KFL1)
      COMMON/OVRINT/IOVEL(2),IOVER(2),IOVEP(2)
      COMMON/FACT/GAM(100)
      COMMON/HOLD/J2STO(3*KFL2-3),J3STO(3*KFL2-3),
     :            J2ANG(3*KFL2),J3ANG(3*KFL2)
      COMMON/INFORM/IREAD,IWRITE,IOUT,ISC0,ISC1,ISC2,ISC3,JSC(3),IALL
      COMMON/DEBUG/IBUG1,IBUG2,IBUG3,NBUG6,NBUG7,IFULL
      COMMON/KRON/IDEL(10,10)
      COMMON/MEDEFN/IHSH,NJ(10),LJ(10),NOSH1(10),NOSH2(10),J1QN1(19,3),
     :     J1QN2(19,3),IJFUL(10)
      COMMON/OVRLAP/MU,NU,MUP,NUP,NONORT,NOVLPS,IROWMU,IROWNU,ICOLMU,
     : ICOLNU,NORTH,IORDER,NCALLS,LMU,LNU,LMUP,LNUP,JMU,JNU,JMUP,JNUP,
     :     IORTH(NORD)
      COMMON/MACHOR/IMATCH(2)
      COMMON/STATES/NCFG,MAXORB,IAJCMP(NWD),LJCOMP(NWD),
     :NJCOMP(NWD),NOCCSH(NCD),NELCSH(5,NCD),NOCORB(5,NCD),J1QNRD(9,NCD)
      COMMON/TERMS/NROWS,ITAB(24),JTAB(24),NTAB(333)
      COMMON/MVALUE/M1,M2,M3,M4,M5,M6,M7,M8,M9,M10,M11,M12,M13,M14,M15,
     :     M16,M17,M18,M19,M20
  300 FORMAT( /15X,
     :   34H COEFFICIENTS OF VARIOUS INTEGRALS//4X,'E-TOTAL(ET)',
     :   '  E-AVERAGE(EAV)','  (ET - EAV) ',
     :   '   INTEGRALS WITH OVERLAPS'/)
  303 FORMAT(7H ISIG =,I3,8H ISIGP =,I3)
  306 FORMAT(5H J1 =,I8,36I3)
  307 FORMAT(3H J2,18X,2HJ3)
  308 FORMAT(3I5,I10,2I5)
  309 FORMAT(4H Y =,F10.6,8H RECUP =,F10.6)
  310 FORMAT(4H Y =,F10.6,8H IDELP =,I3)
  311 FORMAT(3F14.6,4X,1HI,1X,1H(,A3,1H,,A3,1H))
  312 FORMAT(5H Y2 =,F10.6,8H RECUP =,F10.6)
  313 FORMAT(5H Y2 =,F10.6)
  314 FORMAT(3F14.6,4X,1HI,1X,1H(,A3,1H,,A3,1H),1X,
     :     1H<,A3,1H|,A3,3H>**,I2)
  315 FORMAT(3F14.6,4X,1HI,1X,1H(,A3,1H,,A3,1H),1X,
     :     1H<,A3,1H|,A3,3H>**,I2,2H <,A3,1H|,
     :     A3,3H>**,I2)
  316 FORMAT(F14.8,2HL(,A3,I3,1H,,A3,I3,1H):2(1H<,A3,1H|,A3,1H>,I2:))
  317 FORMAT(' IN MATRIX ELEMENT <',I3,'|H0|',I3,'> THERE IS A NON-',
     :  'ZERO COEFFICIENT OF THE'/' ONE-ELECTRON INTEGRAL BETWEEN NON-',
     :  'ORTHOGONAL ORBITALS ',A3,' AND ',A3/' CONSEQUENTLY THE TWO',
     :  ' CONFIGURATIONS ARE NOT ORTHOGONAL AS REQUIRED BY THE PRESENT',
     :  'CODE.')
   10 FORMAT(8H COEFP =,F15.9)
      IZERO = 0
      IONE=1
      IBACK1=0
      IBACK2=0
      IBACK3=0
      IBK = 0
      ISUM=0
      RML=0.D0
      RML2=0.D0
      NOVSTO=NOVLPS
      B1ELC(1) = 0.D0
      B1ELC(2) = 0.D0
      IF (IFULL .NE. 0) WRITE(IWRITE,300)
      IF(IDIAG .NE.1) GO TO 1
*
* --- DIAGONAL HAMILTONIAN MATRIX ELEMENT
*
      Y=0.D0
      Y2=0.D0
      DO 302 J=1,IHSH
      ISIG=IJFUL(J)
      ISIGP=ISIG
      X=NOSH1(J)
      A=X+Y
      IF (DABS(A) .GT. 1.D-10 .AND. IFULL .NE. 0)
     :    WRITE(IWRITE,311) A,X,Y,IAJCMP(ISIG),IAJCMP(ISIGP)
      IF (IALL .NE. 0 .AND. DABS(A) .GT. 1.D-10)
     :   CALL SAVE(4,-A/2.,0,0,ISIG,0,ISIGP,JA,JB,0)
 302  CONTINUE
      IF (NCLOSD .EQ. 0) GO TO 7
      DO 301 J = 1,NCLOSD
      L = LJCLSD(J)
      A = 4*L + 2
      IF (IFULL .NE. 0)
     :  WRITE(IWRITE,311) A, A, Y, IAJCLD(J), IAJCLD(J)
301   CONTINUE
    7 RETURN
*
* --- OFF-DIAGONAL HAMILTONIAN MATRIX ELEMENT
*
    1 CALL ALLADD(IHSH,M3,M4,M5,M6,M7,M8,M9,M10,
     :M11,M12,M13,M14,M15,M16,M17,M18)
*
*     TEST THAT FINAL ANGULAR MOMENTA ARE EQUAL
*
      DO 8 K=2,3
      IF(J1QN1(M6,K).NE.J1QN2(M6,K)) GO TO 7
    8 CONTINUE
*
* --- DETERMINE INTERACTING SHELLS,  ISIG ON L.H.S., ISIGP ON R.H.S.,
*     FOR NON-ZERO 1.D0-ELECTRON MATRIX ELEMENT,  N-1  ELECTRONS MUST BE
*     COMMON TO BOTH SIDES.  THUS THE SUM OF  N(I) = NOSH1(I)-NOSH2(I),
*     I=1,IHSH   MUST BE EQUAL TO  0  OR  2 .   THUS AT NO STAGE CAN
*     N(I) BE GREATER THAN  1 .  IF THIS SUM IS ZERO, THE TWO
*     CONFIGURATIONS ARE MADE UP FROM THE SAME ELECTRONS, WITH TWO
*     DIFFERENT COUPLING SCHEMES.  SINCE THE SPHERICAL HARMONICS ARE
*     EIGENFUNCTIONS OF  DEL**2 ,  THE ORTHOGONALITY OF THE TWO COUPLING
*     SCHEMES WILL BE MAINTAINED AND ORTHOGONALITY GIVES A ZERO RESULT.
*
      IX=0
      IY=0
      DO 9 I=1,IHSH
      IF(I.EQ.MU.OR.I.EQ.NU.OR.I.EQ.MUP.OR.I.EQ.NUP) GO TO 9
      N=NOSH1(I)-NOSH2(I)
      IF(IABS(N).GT.1) GO TO 7
      IF(N) 11,9,12
   11 ISIGP=I
      IY=IY+1
      GO TO 9
   12 ISIG=I
   13 IX=IX+1
    9 CONTINUE
      IF(IX.GT.1.OR.IY.GT.1) GO TO 7
      IF(NOVLPS.NE.0) THEN
         MUSTO=MU
         NUSTO=NU
         MUPSTO=MUP
         NUPSTO=NUP
         NOVSTO=NOVLPS
         LMUSTO=LMU
         LNUSTO=LNU
         LMUPST=LMUP
         LNUPST=LNUP
      ENDIF
      IF(IX.EQ.0.AND.IY.EQ.0) GO TO 3
      IF(NOVLPS.EQ.0) GO TO 24
      IF(IX.EQ.1.AND.IY.EQ.1) THEN
         MGO=MATCH(IZERO,IZERO,IZERO,IZERO)
         IF(MGO.LT.2) GO TO 1000
         GO TO 24
      ENDIF
      IF(IY.EQ.0) GO TO 21
      IBACK1=1
      ISIG=MU
      IF(MU.EQ.NU.OR.NU.EQ.0) IBACK1=0
      MGO=MATCH(IONE,IZERO,IZERO,IZERO)
      IF(MGO.EQ.2) GO TO 24
      GO TO 1000
   25 IBACK1=0
      IBK = 2
      RML=0.D0
      RML2=0.D0
      ISIG=NU
      MGO=MATCH(IZERO,IONE,IZERO,IZERO)
      IF(MGO.LT.2) GO TO 1000
      GO TO 24
   21 IBACK2=1
      ISIGP=MUP
      IF(MUP.EQ.NUP.OR.NUP.EQ.0) IBACK2=0
      MGO=MATCH(IZERO,IZERO,IONE,IZERO)
      IF(MGO.EQ.2) GO TO 24
      GO TO 1000
   22 IBACK2=0
      IBK = 2
      RML=0.D0
      RML2=0.D0
      ISIGP=NUP
      MGO=MATCH(IZERO,IZERO,IZERO,IONE)
      IF(MGO.LT.2) GO TO 1000
   24 IF(NOVLPS.LT.2) THEN
         NONORT=0
      ELSE
         NONORT=1
      ENDIF
      GO TO 4
    3 ISUM=1
      ISIG=1
      ISIGP=1
 1509 IF (NOSH1(ISIG).EQ.0) GO TO 1000
      IF(NOVLPS.EQ.0) GO TO 4
      IF(NOVLPS.EQ.1) GO TO 508
  509 IF(ISIG.NE.MU.AND.ISIG.NE.NU) THEN
         MGO=MATCH(IZERO,IZERO,IZERO,IZERO)
         IF(MGO.EQ.0) THEN
            ISIG=0
            GO TO 1000
         ELSEIF(MGO.EQ.1) THEN
            GO TO 1000
         ELSE
            GO TO 4
         ENDIF
      ELSE
         IF(IBACK3.EQ.0) THEN
            ISIGP=MUP
            IF(MUP.NE.NUP) IBACK3=1
            IF(ISIG.EQ.MU) THEN
               MGO=MATCH(IONE,IZERO,IONE,IZERO)
            ELSE
               MGO=MATCH(IZERO,IONE,IONE,IZERO)
            ENDIF
            IF(MGO.EQ.0) THEN
               ISIG=0
               GO TO 1000
            ELSEIF(MGO.EQ.1) THEN
               IF(IBACK3.EQ.1) THEN
                  MU=MUSTO
                  NU=NUSTO
                  MUP=MUPSTO
                  NUP=NUPSTO
                  NOVLPS=NOVSTO
                  LMU=LMUSTO
                  LNU=LNUSTO
                  LMUP=LMUPST
                  LNUP=LNUPST
                  GO TO 507
               ELSE
                  GO TO 1000
               ENDIF
            ELSE
               GO TO 4
            ENDIF
         ENDIF
  507    IF (IBACK3 .EQ. 1) THEN
            ISIGP=NUP
            IBACK3=0
            RML=0.D0
            RML2=0.D0
            IF(ISIG.EQ.MU) THEN
               MGO=MATCH(IONE,IZERO,IZERO,IONE)
            ELSE
               MGO=MATCH(IZERO,IONE,IZERO,IONE)
            ENDIF
            IF(MGO.EQ.0) THEN
               ISIG=0
               GO TO 1000
            ELSEIF(MGO.EQ.1) THEN
               GO TO 1000
            ELSE
               GO TO 4
            ENDIF
         ENDIF
      ENDIF
  508 IF(ISIG.NE.MU) THEN
         MGO=MATCH(IZERO,IZERO,IZERO,IZERO)
      ELSE
         ISIGP=MUP
         MGO=MATCH(IONE,IZERO,IONE,IZERO)
      ENDIF
      IF(MGO.EQ.0) THEN
         ISIG=0
         GO TO 1000
      ELSEIF(MGO.EQ.1) THEN
         GO TO 1000
      ELSE
         GO TO 4
      ENDIF
    4 IF(IBUG1.GT.0)WRITE(IWRITE,303) ISIG,ISIGP
      LSIG=LJ(ISIG)
      LSIGP=LJ(ISIGP)
*
*     THE ANGULAR MOMENTUM OF THE INTERACTING ELECTRONS MUST BE EQUAL
*
      IF(LSIG.EQ.LSIGP) GO TO 93
      GO TO 1000
*
*     THE SPECTATOR SHELLS MUST HAVE MATCHING QUANTUM NUMBERS
*
   93 DO 16 J=1,IHSH
      IF(J.NE.ISIG) THEN
         DO 253 K=1,3
         JBAR(J,K)=J1QN1(J,K)
  253    CONTINUE
      ENDIF
      IF(J.NE.ISIGP) THEN
         DO 254 K=1,3
         JPBAR(J,K)=J1QN2(J,K)
  254    CONTINUE
      ENDIF
      IF(J.EQ.MU.OR.J.EQ.NU.OR.J.EQ.MUP.OR.J.EQ.NUP) GO TO 16
      IF(J.EQ.ISIG.OR.J.EQ.ISIGP) GO TO 16
      DO 19 K=1,3
      IF(JBAR(J,K).NE.JPBAR(J,K)) GO TO 1000
   19 CONTINUE
   16 CONTINUE
      DO 105 J = 1,IHSH
      NBAR(J)=NOSH1(J)-IDEL(J,ISIG)
  105 CONTINUE
      IF(MUP.NE.0) NBAR(MUP)=NOSH2(MUP)-IDEL(MUP,ISIGP)
      IF(NUP.NE.0) NBAR(NUP)=NOSH2(NUP)-IDEL(NUP,ISIGP)
      IF(NOVLPS.LT.2) GO TO 60
      IF(MU.EQ.NU) THEN
         NKAP=NBAR(MUP)*NBAR(NUP)
      ELSE
         NKAP=NBAR(MU)*NBAR(NU)
      ENDIF
*
* --- INCLUDE MULTIPLICATIVE FACTORS
*
   60 IDELP=0
      IF(ISIG.EQ.IHSH) GO TO 65
      JSIG=ISIG+1
      DO 67 J=JSIG,IHSH
      IF(J.EQ.MUP.OR.J.EQ.NUP) GO TO 67
      IDELP=IDELP+NBAR(J)
   67 CONTINUE
   65 IF(ISIGP.EQ.IHSH) GO TO 170
      JSIGP=ISIGP+1
      DO 68 J=JSIGP,IHSH
      IF(J.EQ.MU.OR.J.EQ.NU) GO TO 68
      IDELP=IDELP+NBAR(J)
   68 CONTINUE
  170 IF(NOVLPS.EQ.1) GO TO 925
      IF (MU.EQ.NU) GO TO 904
      IF (MUP.EQ.NUP) GO TO 905
      IF (NBAR(MU).EQ.NBAR(NUP)) GO TO 906
      MU1N = MIN0(MU,MUP) + 1
      MU1X = MAX0(MU,MUP) - 1
      MU2N = MIN0(NU,NUP) + 1
      MU2X = MAX0(NU,NUP) - 1
      MU1=MU
      NU1=NU
      GO TO 933
  906 MU1N = MIN0(MU,NUP) + 1
      MU1X = MAX0(MU,NUP) - 1
      MU2N = MIN0(NU,MUP) + 1
      MU2X = MAX0(NU,MUP) - 1
      MU1=MU
      NU1=NU
      GO TO 933
  905 MU1N = MIN0(MU,MUP) + 1
      MU1X = MAX0(MU,MUP) - 1
      MU2N = MIN0(MUP,NU) + 1
      MU2X = MAX0(MUP,NU) - 1
      MU1=MU
      NU1=NU
      GO TO 933
  904 MU1N = MIN0(MU,MUP) + 1
      MU1X = MAX0(MU,MUP) - 1
      MU2N = MIN0(MU,NUP) + 1
      MU2X = MAX0(MU,NUP) - 1
      MU1=MUP
      NU1=NUP
  933 IF (MU1N .GT. MU1X) GO TO 934
      DO 962 J = MU1N,MU1X
      IF (J.EQ.MU .OR. J.EQ.NU .OR. J.EQ.MUP .OR. J.EQ.NUP) GO TO 962
      IDELP = IDELP + NBAR(J)*NBAR(MU1)
  962 CONTINUE
  934 IF (MU2N .GT. MU2X) GO TO 70
      DO 963 J = MU2N,MU2X
      IF (J.EQ.MU .OR. J.EQ.NU .OR. J.EQ.MUP .OR. J.EQ.NUP) GO TO 963
      IDELP = IDELP + NBAR(J)*NBAR(NU1)
  963 CONTINUE
      GO TO 70
  925 MUMIN1=MIN0(MU,MUP)+1
      MUMAX1=MAX0(MU,MUP)-1
      IF(MUMIN1.GT.MUMAX1) GO TO 70
      DO 927 J=MUMIN1,MUMAX1
      IDELP=IDELP+NBAR(J)*NBAR(MU)
  927 CONTINUE
   70 FACTR=(-1.D0)**IDELP*SQRT(DFLOAT(NOSH1(ISIG)*NOSH2(ISIGP)))
      IF(NOVLPS.LT.2) GO TO 102
      IF(MU.NE.NU.AND.MUP.NE.NUP) GO TO 102
      N1=IOVEP(1)
      N2=IOVEP(2)
      N3=N1+N2
      FACTR=FACTR*EXP((GAM(N3+1)-GAM(N1+1)-GAM(N2+1))/2.D0)
*
* --- SET UP  J2  AND  J3  ARRAYS
*
  102 M1=IHSH-2
      M2=M6-2
      J2(1,1)=ISIG
      J2(1,2)=M11
      J2(1,3)=M9
      J3(1,1)=ISIGP
      J3(1,2)=M11
      J3(1,3)=M10
      IF(ISIG.EQ.1) GO TO 29
      J2(2,1)=1
      GO TO 30
   29 J2(2,1)=M9
   30 IF(ISIG.EQ.2) GO TO 32
      J2(2,2)=2
      GO TO 33
   32 J2(2,2)=M9
   33 J2(2,3)=M4
      IF(ISIGP.EQ.1) GO TO 35
      J3(2,1)=1
      GO TO 36
   35 J3(2,1)=M10
   36 IF(ISIGP.EQ.2) GO TO 38
      J3(2,2)=2
      GO TO 39
   38 J3(2,2)=M10
   39 J3(2,3)=M7
      IF(IHSH.LT.3) GO TO 40
      DO 42 J=3,IHSH
      J2(J,1)=M1+J
      J2(J,3)=M1+J+1
      J3(J,1)=M2+J
      IF(J.EQ.IHSH) GO TO 44
      J3(J,3)=M2+J+1
      GO TO 45
   44 J3(J,3)=M1+J+1
   45 IF(J.EQ.ISIG) GO TO 47
      J2(J,2)=J
      GO TO 48
   47 J2(J,2)=M9
   48 IF(J.EQ.ISIGP) GO TO 50
      J3(J,2)=J
      GO TO 42
   50 J3(J,2)=M10
   42 CONTINUE
*
* --- STORE  J2  AND  J3  ARRAYS FOR USE IN SPIN RECOUPLING COEFFICIENT
*
   40 I1=2
      NJ1S=M11
      NJ23S=M4
*
* --- Added for NJGRAF
*
      DO 100 j = 1,M11
         FREE(j) = .FALSE.
  100 CONTINUE
*
   17 IF(NOVLPS.NE.0) CALL CNDENS(I1)
      I1=0
      NJ23=NJ23S-1
      DO 51 J=1,NJ23
      DO 52 K=1,3
      I1=I1+1
      J2STO(I1)=J2(J,K)
      J3STO(I1)=J3(J,K)
   52 CONTINUE
   51 CONTINUE
*
*   SUM OVER PARENTS OF ISIG AND ISIGP
*
      N1=NOSH1(ISIG)
      L1=LSIG+1
      K1=NTAB1(N1,L1)
      KK1=ITAB(K1)
      N2=NOSH2(ISIGP)
      L2=LSIGP+1
      K2=NTAB1(N2,L2)
      KK2=ITAB(K2)
      DO 111 JJ1=1,KK1
*
*    CHECK ON TRIANGULAR CONDITIONS
*
      IN3=2*LSIG
      IJK1=3*(JJ1-1)+JTAB(K1)
      DO 113 KK=2,3
      IN1=NTAB(IJK1+KK)-1
      IN2=J1QN1(ISIG,KK)-1
      IF(IN1.GT.(IN2+IN3)) GO TO 111
      IF(IN1.LT.IABS(IN2-IN3)) GO TO 111
      IN3=1
  113 CONTINUE
      DO 112 JJ2=1,KK2
      IN3=2*LSIGP
      IJK2=3*(JJ2-1)+JTAB(K2)
      DO 114 KK=2,3
      IN1=NTAB(IJK2+KK)-1
      IN2=J1QN2(ISIGP,KK)-1
      IF(IN1.GT.(IN2+IN3)) GO TO 112
      IF(IN1.LT.IABS(IN2-IN3)) GO TO 112
      IN3=1
  114 CONTINUE
      DO 115 KK=1,3
      JBAR(ISIG,KK)=NTAB(IJK1+KK)
      JPBAR(ISIGP,KK)=NTAB(IJK2+KK)
  115 CONTINUE
      IF(ISIG.EQ.MU.OR.ISIG.EQ.NU) GO TO 116
      DO 117 KK=1,3
      IF(JBAR(ISIG,KK).NE.JPBAR(ISIG,KK)) GO TO 112
  117 CONTINUE
  116 IF(ISIGP.EQ.MUP.OR.ISIGP.EQ.NUP) GO TO 118
      DO 119 KK=1,3
      IF(JBAR(ISIGP,KK).NE.JPBAR(ISIGP,KK)) GO TO 112
  119 CONTINUE
  118 IF(NOVLPS.EQ.1) THEN
         DO 120 KK=1,3
         IF(JBAR(MU,KK).NE.JPBAR(MUP,KK)) GO TO 112
  120    CONTINUE
      ENDIF
*
* --- CALCULATE FRACTIONAL PARENTAGE COEFFICIENTS
*
      Y=1.D0
      Y2=0.D0
      IF(LSIG.EQ.0) GO TO 26
      N=NOSH1(ISIG)
      IVI=J1QN1(ISIG,1)
      ILI=(J1QN1(ISIG,2)-1)/2
      ISI=J1QN1(ISIG,3)
      IVJ=JBAR(ISIG,1)
      ILJ=(JBAR(ISIG,2)-1)/2
      ISJ=JBAR(ISIG,3)
      CALL CFP(LSIG,N,IVI,ILI,ISI,IVJ,ILJ,ISJ,COEFP)
      IF(IBUG1.GT.0) WRITE(IWRITE,10) COEFP
      Y=Y*COEFP
      N=NOSH2(ISIGP)
      IVI=J1QN2(ISIGP,1)
      ILI=(J1QN2(ISIGP,2)-1)/2
      ISI=J1QN2(ISIGP,3)
      IVJ=JPBAR(ISIGP,1)
      ILJ=(JPBAR(ISIGP,2)-1)/2
      ISJ=JPBAR(ISIGP,3)
      CALL CFP(LSIG,N,IVI,ILI,ISI,IVJ,ILJ,ISJ,COEFP)
      IF(IBUG1.GT.0) WRITE(IWRITE,10) COEFP
      Y=Y*COEFP
   26 IF(NOVLPS.EQ.2) Y2=Y
      IF(DABS(Y).LT.1.D-14.AND.DABS(Y2).LT.1.D-14) GO TO 112
*
* --- ORBITAL RECOUPLING COEFFICIENT
*
      J1(M11)=LSIG+LSIG+1
      K=2
*
* --- SET  J1  ARRAY
*
   64 DO 53 J=1,IHSH
      J1(J)=JBAR(J,K)
      IF (J.EQ.MUP .OR. J.EQ.NUP) J1(J) = JPBAR(J,K)
   53 CONTINUE
      DO 56 J=M4,M6
      J1(J)=J1QN1(J,K)
   56 CONTINUE
      DO 57 J=M7,M8
      J1(J)=J1QN2(J-M3,K)
   57 CONTINUE
      J1(M9)=J1QN1(ISIG,K)
      J1(M10)=J1QN2(ISIGP,K)
*
      IF(IBUG1.LT.1.OR.IBUG3.EQ.1) GO TO 77
      WRITE(IWRITE,306) (J1(J),J=1,M11)
      WRITE(IWRITE,307)
      DO 80 J=1,NJ23
      WRITE(IWRITE,308) (J2(J,KL),KL=1,3),(J3(J,KL),KL=1,3)
   80 CONTINUE
*
* --- CALCULATE RECOUPLING COEFFICIENT
*
   77 IF(NOVLPS.LT.2.OR.MU.EQ.NU.OR.MUP.EQ.NUP) GO TO 78
      IF(IMATCH(1).EQ.0) GO TO 79
      IF(LMU.NE.LMUP.OR.LNU.NE.LNUP) GO TO 79
      DO 937 KK=1,3
      IF(JBAR(MU,KK).NE.JPBAR(MUP,KK)) GO TO 79
      IF(JBAR(NU,KK).NE.JPBAR(NUP,KK)) GO TO 79
  937 CONTINUE
     
   78 CALL NJGRAF(RECUP,FAIL)
      GO TO 81
   79 RECUP=0.D0
   81 Y=Y*RECUP
      IF(IBUG1.GT.0) WRITE(IWRITE,309) Y,RECUP
      IF(NOVLPS.LT.2) GO TO 23
      IF(IMATCH(2).EQ.0) GO TO 84
      IF(MU.EQ.NU.OR.MUP.EQ.NUP) GO TO 84
      IF(LMU.NE.LNUP.OR.LNU.NE.LMUP) GO TO 84
      DO 938 KK=1,3
      IF(JBAR(MU,KK).NE.JPBAR(NUP,KK)) GO TO 84
      IF(JBAR(NU,KK).NE.JPBAR(MUP,KK)) GO TO 84
  938 CONTINUE
      I1=0
      DO 82 J=1,NJ23
      DO 83 KK=1,3
      I1=I1+1
      J2(J,KK)=J2STO(I1)
      J3(J,KK)=J3STO(I1)
   83 CONTINUE
   82 CONTINUE
      JSTO=J3(IROWMU,ICOLMU)
      J3(IROWMU,ICOLMU)=J3(IROWNU,ICOLNU)
      J3(IROWNU,ICOLNU)=JSTO
     
      CALL NJGRAF(RECUP,FAIL)
   86 Y2=Y2*RECUP
   89 IF(IBUG1.GT.0) WRITE(IWRITE,312) Y2,RECUP
      GO TO 23
   84 RECUP=0.D0
      GO TO 86
*
* --- SPIN RECOUPLING COEFFICIENT
*
  23  I1=0
      DO 62 J=1,NJ23
      DO 63 KK=1,3
      I1=I1+1
      J2(J,KK)=J2STO(I1)
      J3(J,KK)=J3STO(I1)
   63 CONTINUE
   62 CONTINUE
      IF(K.EQ.3) GO TO 101
      J1(M11 )=2
      K=3
      GO TO 64
  101 VALML=Y
      IF(NOVLPS.EQ.2) VALML2=Y2
      RML=RML+VALML
      IF(NOVLPS.EQ.2) RML2=RML2+VALML2
  112 CONTINUE
  111 CONTINUE
      Y=RML*FACTR
      IF(IBUG1.GT.1) WRITE(IWRITE,310) Y,IDELP
      IF(NOVLPS.LT.2) GO TO 90
      Y2=RML2*FACTR*(-1.D0)**NKAP
      IF(IBUG1.GT.1) WRITE(IWRITE,313) Y2
   90 X=0.D0
      A=X+Y
      JSIG=IJFUL(ISIG)
      JSIGP=IJFUL(ISIGP)
      JSL = MIN0(JSIG,JSIGP)
      JSU = MAX0(JSIG,JSIGP)
      JS1 = ((JSU-1)*(JSU-2))/2 + JSL
      JS2 = IORTH(JS1)
      IF ( (DABS(Y).GT.1.D-10 .OR. DABS(Y2).GT.1.D-10)  .AND.
     :     (JS2.EQ.1 .OR. JSIG.EQ.JSIGP) )  THEN
           WRITE(IWRITE,317) JA,JB,IAJCMP(JSIG),IAJCMP(JSIGP)
           STOP
      END IF
  110 IF (NOVLPS.EQ.0) GO TO 98
      JMU = IJFUL(MU)
      JMUP = IJFUL(MUP)
   96 IF (NOVLPS .EQ. 1) GO TO 91
      JNU = IJFUL(NU)
      JNUP = IJFUL(NUP)
      GO TO 92
   98 IF (DABS(Y) .GT. 1.D-10 ) THEN
         IF (IFULL .NE. 0)
     :       WRITE(IWRITE,311) A,X,Y,IAJCMP(JSIG),IAJCMP(JSIGP)
         B1ELC(IBK+1) = Y
         YY = -Y
         IF (IALL.NE.0) YY = YY/2.
         CALL SAVE(4,YY,0,0,JSIG,0,JSIGP,JA,JB,0)
       ENDIF
       GO TO 1000
   91  IF (DABS(Y) .GE. 1.D-10 ) THEN
         IF (IFULL .NE. 0)
     :      WRITE(IWRITE,314) A,X,Y,IAJCMP(JSIG),
     :     IAJCMP(JSIGP),IAJCMP(JMU),IAJCMP(JMUP),IOVEP(1)
         YY = -Y
         IF (IALL.NE.0) YY = YY/2.
         CALL OVLSET(1,JMU,JMUP,IOVEP(1),0,0,0,IPTR)
         CALL SAVE(4,YY,0,0,JSIG,0,JSIGP,JA,JB,IPTR)
         B1ELC(IBK+1) = Y
      ENDIF
      GO TO 1000
   92 IF (DABS(Y) .GE. 1.D-10) THEN
         IF (IFULL .NE. 0)
     :      WRITE(IWRITE,315) A,X,Y,IAJCMP(JSIG),
     :      IAJCMP(JSIGP),IAJCMP(JMU),IAJCMP(JMUP),IOVEP(1),IAJCMP(JNU),
     :      IAJCMP(JNUP),IOVEP(2)
         YY = -Y
         IF (IALL.NE.0) YY = YY/2.
         CALL OVLSET(2,JMU,JMUP,IOVEP(1),JNU,JNUP,IOVEP(2),IPTR)
         CALL SAVE(4,YY,0,0,JSIG,0,JSIGP,JA,JB,IPTR)
         B1ELC(IBK+1) = Y
      ENDIF
      IF (DABS(Y2) .LT. 1.D-10) GO TO 1000
      A=X+Y2
       IF (IFULL .NE. 0) WRITE(IWRITE,315) A,X,Y2,IAJCMP(JSIG),
     :   IAJCMP(JSIGP),IAJCMP(JNU),IAJCMP(JMUP),IOVEP(2),IAJCMP(JMU),
     :   IAJCMP(JNUP),IOVEP(1)
         YY = -Y2
         IF (IALL.NE.0) YY = YY/2.
         CALL OVLSET(2,JMU,JNUP,IOVEP(1),JNU,JMUP,IOVEP(2),IPTR)
         CALL SAVE(4,YY,0,0,JSIG,0,JSIGP,JA,JB,IPTR)
         B1ELC(IBK+2) = Y2
 1000 IF(NOVSTO.EQ.0) GO TO 289
      MU=MUSTO
      NU=NUSTO
      MUP=MUPSTO
      NUP=NUPSTO
      NOVLPS=NOVSTO
      LMU=LMUSTO
      LNU=LNUSTO
      LMUP=LMUPST
      LNUP=LNUPST
      IF(IBACK1.EQ.1) GO TO 25
      IF(IBACK2.EQ.1) GO TO 22
  289 RML=0.D0
      RML2=0.D0
      IF(ISUM.EQ.0) RETURN
      IF (MGO .EQ. 0) RETURN
      IF(IBACK3.EQ.1) GO TO 509
      ISIG=ISIG+1
      ISIGP=ISIG
      IF(ISIG.LE.IHSH) GO TO 1509
      RETURN
      END
*
*     ------------------------------------------------------------------
*       J 2 3 A N G
*     ------------------------------------------------------------------
*
      SUBROUTINE J23ANG(IRHO,ISIG,IRHOP,ISIGP,JANGDI)
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      PARAMETER (NWD=30,NCD=100,NORD=NWD*(NWD-1)/2)
*
      COMMON/INFORM/IREAD,IWRITE,IOUT,ISC0,ISC1,ISC2,ISC3,JSC(3),IALL
      COMMON/DEBUG/IBUG1,IBUG2,IBUG3,NBUG6,NBUG7,IFULL
      COMMON/MVALUE/M1,M2,M3,M4,M5,M6,M7,M8,M9,M10,M11,M12,M13,M14,M15,
     :     M16,M17,M18,M19,M20
      PARAMETER(KFL1=60,KFL2=12)
      LOGICAL FREE
      COMMON/COUPLE/NJ1S,NJ23S,J1(KFL1),J2(KFL2,3),J3(KFL2,3),FREE(KFL1)
      COMMON/HOLD/J2SPIN(3*KFL2-3),J3SPIN(3*KFL2-3),
     :            J2ANG(3*KFL2),J3ANG(3*KFL2)
      COMMON/OVRLAP/MU,NU,MUP,NUP,NONORT,NOVLPS,IROWMU,IROWNU,ICOLMU,
     : ICOLNU,NORTH,IORDER,NCALLS,LMU,LNU,LMUP,LNUP,JMU,JNU,JMUP,JNUP,
     :     IORTH(NORD)
      COMMON/J23NOR/NJ23SP,NJ23AN
*
* === SETS UP J2 AND J3 ARRAYS FOR DIRECT ANGULAR INTEGRAL CALL OF NJSYM
*
  303 FORMAT(3H J2,18X,2HJ3)
  304 FORMAT(3I5,I10,2I5)
*
*     HAVE THE J2 AND J3 ARRAYS ALREADY BEEN  SET.  IF NOT, THEN GO TO 2
*
      IF(JANGDI) 2,2,1
*
* --- ROWS 3 TO M4 OF SPIN J2 AND J3 ARE SAME AS ROWS 4 TO (M4+1) OF
*     ANGULAR J2 AND J3
*
    2 I1=6
*
*     CURRENT VALUE OF NJ23S IS FOR RECOUPLING COEFFICIENT
*     AS SET IN J1
*
      NJ23=NJ23S-2
      DO 103 J=3,NJ23
      JP1=J+1
      DO 104 K=1,3
      I1=I1+1
      J2(JP1,K)=J2SPIN(I1)
      J3(JP1,K)=J3SPIN(I1)
  104 CONTINUE
  103 CONTINUE
*
* --- SET ROWS 1, 2 AND 3
*
      IF(M1) 105,106,105
  105 J2(3,1)=ISIG
      GO TO 107
  106 J2(3,1)=M9
  107 IF(M2) 109,110,109
  109 IF(ISIGP.EQ.MUP .AND. MU.NE.NU) GO TO 11
      IF(ISIGP.EQ.NUP .AND. MU.NE.NU) GO TO 12
      J3(3,1)=ISIGP
      GO TO 111
   11 J3(3,1)=MU
      GO TO 111
   12 J3(3,1)=NU
      GO TO 111
  110 J3(3,1)=M11
  111 J2(2,3)=M9
      J2(2,1)=IRHO
      J2(2,2)=M13
      J2(1,3)=M14
      J2(3,2)=M14
      J2(3,3)=M10
      J2(1,1)=M16
      J2(1,2)=M17
      J3(3,2)=M16
      J3(3,3)=M12
      J3(1,2)=M13
      J3(1,1)=M17
      J3(1,3)=M15
      J3(2,3)= M11
      IF(IRHOP.EQ.MUP .AND. MU.NE.NU) GO TO 21
      IF(IRHOP.EQ.NUP .AND. MU.NE.NU) GO TO 22
      J3(2,1)=IRHOP
      GO TO 23
   21 J3(2,1)=MU
      GO TO 23
   22 J3(2,1)=NU
   23 J3(2,2)=M15
*
* --- STORE J2 AND J3 FOR USE IN CALCULATING THE EXCHANGE TERM
*
      NJ23=NJ23S-1
      IF (NOVLPS .LT. 2) GO TO 6
      IF (MU.EQ.NU .AND. MUP.EQ.NUP) GO TO 6
      IF (MU.NE.NU .AND. MUP.NE.NUP) GO TO 6
      IF (MU .NE.NU .AND. MUP.EQ.NUP) GO TO 7
      DO 8 J = 1,3
      I2 = 4-J
      I1 = I2+1
      DO 9 K = 1,3
      J2(I1,K) = J2(I2,K)
    9 CONTINUE
    8 CONTINUE
      J2(1,1) = MUP
      J2(1,2) = NUP
      J2(1,3) = MU
      GO TO 6
    7 DO 15 J = 1,3
      I2 = 4 - J
      I1 = I2+1
      DO 16 K =1,3
      J3(I1,K) = J3(I2,K)
   16 CONTINUE
   15 CONTINUE
      J3(1,1) = MU
      J3(1,2) = NU
      J3( 1,3) = MUP
    6 I1=0
      DO 535 J=1,NJ23
      DO 536 K=1,3
      I1=I1+1
      J2ANG(I1)=J2(J,K)
      J3ANG(I1)=J3(J,K)
  536 CONTINUE
  535 CONTINUE
      NJ23AN=NJ23+1
      JANGDI=1
    3 IF(IBUG2-1) 209,206,206
*
*     PRINT-OUT OF VALUES IN NJSYM  IF IBUG3=1
*
  206 IF(IBUG3-1) 207,209,207
  207 WRITE(IWRITE,303)
      DO 208 J=1,NJ23
      WRITE(IWRITE,304) (J2(J,K),K=1,3),(J3(J,K),K=1,3)
  208 CONTINUE
  209 RETURN
*
* --- SET J2 AND J3 ARRAYS FROM STORE OF PREVIOUS CALCULATIONS
*
    1 I1=0
      NJ23S=NJ23AN
      NJ23=NJ23S-1
      DO 4 J=1,NJ23
      DO 5 K=1,3
      I1=I1+1
      J2(J,K)=J2ANG(I1)
      J3(J,K)=J3ANG(I1)
    5 CONTINUE
    4 CONTINUE
      GO TO 3
      END
*
*     ------------------------------------------------------------------
*       J 2 3 S P N
*     ------------------------------------------------------------------
*
      SUBROUTINE J23SPN(IRHO,ISIG,IRHOP,ISIGP,JSNDIR)
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      PARAMETER (NWD=30,NCD=100,NORD=NWD*(NWD-1)/2)
*
      COMMON/INFORM/IREAD,IWRITE,IOUT,ISC0,ISC1,ISC2,ISC3,JSC(3),IALL
      COMMON/DEBUG/IBUG1,IBUG2,IBUG3,NBUG6,NBUG7,IFULL
      COMMON/MVALUE/M1,M2,M3,M4,M5,M6,M7,M8,M9,M10,M11,M12,M13,M14,M15,
     :     M16,M17,M18,M19,M20
      PARAMETER(KFL1=60,KFL2=12)
      LOGICAL FREE
      COMMON/COUPLE/NJ1S,NJ23S,J1(KFL1),J2(KFL2,3),J3(KFL2,3),FREE(KFL1)
      COMMON/HOLD/J2SPIN(3*KFL2-3),J3SPIN(3*KFL2-3),
     :            J2ANG(3*KFL2),J3ANG(3*KFL2)
      COMMON/OVRLAP/MU,NU,MUP,NUP,NONORT,NOVLPS,IROWMU,IROWNU,ICOLMU,
     : ICOLNU,NORTH,IORDER,NCALLS,LMU,LNU,LMUP,LNUP,JMU,JNU,JMUP,JNUP,
     :     IORTH(NORD)
      COMMON/J23NOR/NJ23SP,NJ23AN
*
* === SET UP THE J2 AND J3 ARRAYS FOR THE DIRECT SPIN INTEGRAL CALL
*     OF NJSYM
*
  303 FORMAT(3H J2,18X,2HJ3)
  304 FORMAT(3I5,I10,2I5)
*
*     HAVE THE J2 AND J3 ARRAYS ALREADY BEEN  SET.  IF NOT, THEN GO TO 2
*
      IF(JSNDIR) 2,2,1
*
* --- SET THIRD ROW OF J2 AND J3
*
    2 IF(IRHO-1) 271,272,271
  271 J2(3,1)=1
      GO TO 273
  272 IF(M1) 274,275,274
  275 J2(3,1)=M10
      GO TO 276
  274 J2(3,1)=M9
      GO TO 276
  273 IF(IRHO-2) 277,278,277
  277 J2(3,2)=2
      GO TO 284
  278 IF(M1) 280,281,280
  280 J2(3,2)=M9
      GO TO 284
  281 J2(3,2)=M10
      GO TO 284
  276 IF(ISIG-2) 277,281,277
  284 J2(3,3)=M4
      IF(IRHOP-1) 285,286,285
  285 J3(3,1)=1
      GO TO 287
  286 IF(M2) 288,289,288
  288 J3(3,1)=M11
      GO TO 290
  289 J3(3,1)=M12
      GO TO 290
  287 IF(IRHOP-2) 291,292,291
  291 J3(3,2)=2
      GO TO 293
  292 IF(M2) 294,295,294
  295 J3(3,2)=M12
      GO TO 293
  294 J3(3,2)=M11
      GO TO 293
  290 IF(ISIGP-2) 291,295,291
  293 J3(3,3)=M7
*
* --- SET ROWS 4,5,.. ETC.
*
      IF(M4-4) 203,202,202
  202 DO 470 J=4,M4
      J2(J,1)=M4 +J-4
      J2(J,3)=M4+J-3
      IF(ISIG+1-J) 471,472,471
  471 IF(M1) 473,474,473
  473 IF(IRHO+1-J) 474,475,474
  474 J2(J,2)=J-1
      GO TO 476
  472 J2(J,2)=M10
      GO TO 476
  475 J2(J,2)=M9
  476 J3(J,1)=M7+J-4
      IF(J-M4 ) 482,483,482
  483 J3(J,3)=J2(J,3)
      GO TO 484
  482 J3(J,3)=M7+J-3
  484 IF(ISIGP+1-J) 477,478,477
  477 IF(M2) 479,480,479
  479 IF(IRHOP+1-J) 480,481,480
  480 J3(J,2)=J-1
      GO TO 470
  478 J3(J,2)=M12
      GO TO 470
  481 J3(J,2)=M11
  470 CONTINUE
*
* --- SET FIRST TWO ROWS, CORRESPONDING TO COUPLING OF INTERACTING
*     ELECTRONS WITHIN THEIR SHELLS
*
  203 J2(2,3)=M10
      J2(1,2) = M13
      J2(2,2) = M14
      J2(1,3) = M9
      IF(M1) 82,83,82
   82 J2(1,1) = IRHO
      J2(2,1) = ISIG
      GO TO 84
   83 J2(1,1) = ISIG
      J2(2,1) = M9
   84 J3(2,3) = M12
      J3(1,2) = M13
      J3(2,2) = M14
      J3(1,3) = M11
      IF(M2) 85,86,85
   85 J3(1,1) = IRHOP
      J3(2,1) = ISIGP
      GO TO 187
   86 J3(1,1) = ISIGP
      J3(2,1) = M11
*
* --- STORE J2,J3 ARRAYS FOR USE IN CALCULATING EXCHANGE INTEGRAL
*
  187 I1=2
      IF(NOVLPS.NE.0) CALL CNDENS(I1)
      I1=0
      NJ23=NJ23S-1
      NJ23SP=NJ23S
      DO 451 J=1,NJ23
      DO 452 K=1,3
      I1=I1+1
      J2SPIN(I1)=J2(J,K)
      J3SPIN(I1)=J3(J,K)
  452 CONTINUE
  451 CONTINUE
      JSNDIR=1
    3 IF(IBUG2-1) 570,6,6
*
*     PRINT-OUT OF VALUES IN NJSYM  IF IBUG3=1
*
    6 IF(IBUG3-1) 200,570,200
  200 WRITE(IWRITE,303)
      DO 201 J=1,NJ23
      WRITE(IWRITE,304) (J2(J,K),K=1,3),(J3(J,K),K=1,3)
  201 CONTINUE
  570 RETURN
*
* --- SET J2 AND J3 ARRAYS FROM STORE OF PREVIOUS CALCULATIONS
*
    1 I1=0
      NJ23S=NJ23SP
      NJ23=NJ23S-1
      DO 4 J=1,NJ23
      DO 5 K=1,3
      I1=I1+1
      J2(J,K)=J2SPIN(I1)
      J3(J,K)=J3SPIN(I1)
    5 CONTINUE
    4 CONTINUE
      GO TO 3
      END
*     MCHF_NONH  (Part 2 of 2)
*     ------------------------------------------------------------------
*       M A T C H
*     ------------------------------------------------------------------
*
      INTEGER FUNCTION MATCH(IA, IB, IC, ID)
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      PARAMETER (NWD=30,NCD=100,NORD=NWD*(NWD-1)/2)
*
      COMMON/OVRINT/IOVEL(2),IOVER(2),IOVEP(2)
      COMMON/OVRLAP/MU,NU,MUP,NUP,NONORT,NOVLPS,IROWMU,IROWNU,ICOLMU,
     : ICOLNU,NORTH,IORDER,NCALLS,LMU,LNU,LMUP,LNUP,JMU,JNU,JMUP,JNUP,
     :     IORTH(NORD)
      COMMON/MEDEFN/IHSH,NJ(10),LJ(10),NOSH1(10),NOSH2(10),J1QN1(19,3),
     :     J1QN2(19,3),IJFUL(10)
      COMMON/INFORM/IREAD,IWRITE,IOUT,ISC0,ISC1,ISC2,ISC3,JSC(3),IALL
      COMMON/DEBUG/IBUG1,IBUG2,IBUG3,NBUG6,NBUG7,IFULL
      COMMON/DIAGNL/IDIAG,JA,JB
      COMMON/MACHOR/IMATCH(2)
*
* --- Tries to match the angular momenta of the spectator
*     non-orthogonal orbitals.
*
      IMATCH(1)=0
      IMATCH(2)=0
      IM = 0
      NA = 0
      NB = 0
      NC = 0
      ND = 0
  200 IF (MU .NE. 0) NA = NOSH1(MU)
      IF (NU .NE. 0 .AND. MU .NE. NU) NB = NOSH1(NU)
      IF (MUP .NE. 0) NC = NOSH2(MUP)
      IF (NUP .NE. 0 .AND. MUP .NE. NUP) ND = NOSH2(NUP)
      NA = NA - IA
      NB = NB - IB
      NC = NC - IC
      ND = ND - ID
      IF (NA.LT.0 .OR. NB.LT.0 .OR. NC.LT.0 .OR. ND.LT.0) GO TO 20
      IF ((NA+NB) .EQ. (NC+ND)) GO TO 14
   20 MATCH = 1
      RETURN
   14 IF (NA.EQ.0 .AND. NB.EQ.0 .AND. NC.EQ.0 .AND. ND.EQ.0) GO TO 17
      IF ((NA.EQ.0 .OR. NB.EQ.0).AND.(NC.EQ.0 .OR. ND.EQ.0)) GO TO 15
      NOVLPS=2
      GO TO 16
   15 IF (NA .EQ. 0) MU = NU
      IF (NC .EQ. 0) MUP = NUP
      NU = 0
      NUP = 0
      NOVLPS = 1
      LMU = LJ(MU)
      LMUP = LJ(MUP)
      IF (LMU .NE. LMUP) GO TO 20
      IOVEL(1) = MU
      IOVER(1) = MUP
      IOVEP(1) = NA + NB
      MATCH = 2
      GO TO 30
   17 MU = 0
      NU = 0
      MUP = 0
      NUP = 0
      NOVLPS = 0
      MATCH = 2
      RETURN
   16 IF(NA.EQ.0) MU=NU
      IF(NB.EQ.0) NU=MU
      IF(NC.EQ.0) MUP=NUP
      IF(ND.EQ.0) NUP=MUP
      LA = LJ(MU)
      LB = LJ(NU)
      LC = LJ(MUP)
      LD = LJ(NUP)
      IF (LA .EQ. LB) GO TO 21
      IF (LA .EQ. LC) GO TO 22
      IF (LA .NE. LD .OR. NA .NE. ND) GO TO 20
      IF (LB .NE. LC) GO TO 20
      IOVEL(1) = MU
      IOVER(1) = NUP
      IOVEP(1) = NA
      IOVEL(2) = NU
      IOVER(2) = MUP
      IOVEP(2) = NB
      MATCH = 2
      IM = 1
      GO TO 30
   22 IF (LB.NE.LD .OR. NB.NE.ND) GO TO 20
      IOVEL(1) = MU
      IOVER(1) = MUP
      IOVEP(1) = NA
      IOVEL(2) = NU
      IOVER(2) = NUP
      IOVEP(2) = NB
      MATCH = 2
      GO TO 30
   21 IF (LA.NE.LC .OR. LA.NE.LD) GO TO 20
      IF(MU .EQ. NU) GO TO 23
      IF (MUP .EQ. NUP) GO TO 24
      IF (NA.GT.1 .OR. NB.GT.1 .OR. NC.GT.1 .OR. ND.GT.1) GO TO 25
      IOVEL(1) = MU
      IOVER(1) = MUP
      IOVEP(1) = 1
      IOVEL(2) = NU
      IOVER(2) = NUP
      IOVEP(2) = 1
      MATCH = 2
      GO TO 30
   25 WRITE(IWRITE,300) MU,NU,MUP,NUP,JA,JB
  300 FORMAT('THE FOLLOWING SUBSHELLS HAVE A COMMON L-VALUE BUT',
     :  ' CONTAIN TOO MANY ELECTRONS FOR THIS CODE',5X,4I3/
     :    5X,'THE MATRIX ELEMENT IS (',I2,1H,,I2,1H))
      MATCH = 0
      RETURN
   23 IF ((NC+ND) .GT. 2) GO TO 25
      IOVEL(1) = MU
      IOVER(1) = MUP
      IOVEP(1) = NC
      IOVEL(2) = MU
      IOVER(2) = NUP
      IOVEP(2) = ND
      MATCH = 2
      GO TO 30
   24 IF ((NA+NB) .GT. 2) GO TO 25
      IOVEL(1) = MU
      IOVER(1) = MUP
      IOVEP(1) = NA
      IOVEL(2) = NU
      IOVER(2) = MUP
      IOVEP(2) = NB
      MATCH = 2
*
*      FINAL CHECK ON ALLOWED NON-ORTHOGONALITY
*
   30 K1=1
      K2=1
      IGO=1
   31 I1=IOVEL(K1)
      I2=IOVER(K2)
      I5=IJFUL(I1)
      I6=IJFUL(I2)
      I3=MIN0(I5,I6)
      I4=MAX0(I5,I6)
      I1=(I4-2)*(I4-1)/2+I3
      I2=IORTH(I1)
      GO TO (32,33,35,36),IGO
   32 IF(I2.NE.1) GO TO 34
      IF(NOVLPS.EQ.1) GO TO 33
      K1=2
      K2=2
      IGO=2
      GO TO 31
   33 IF(I2.EQ.1) IMATCH(IM+1)=1
   34 IF(NOVLPS.EQ.1) THEN
         IF(IMATCH(1).EQ.0) GO TO 20
         RETURN
      ENDIF
      K1=1
      K2=2
      IGO=3
      GO TO 31
   35 IF(I2.NE.1) GO TO 37
      K1=2
      K2=1
      IGO=4
      GO TO 31
   36 IF(I2.EQ.1) IMATCH(2-IM)=1
   37 IF(IMATCH(1).EQ.0.AND.IMATCH(2).EQ.0) GO TO 20
      RETURN
      END
*
*     ------------------------------------------------------------------
*       M E K E E P
*     ------------------------------------------------------------------
*
      SUBROUTINE MEKEEP(IRHO,ISIG,IRHOP,ISIGP)
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      PARAMETER (NWD=30,NCD=100,NORD=NWD*(NWD-1)/2)
*
      COMMON/MEDEFN/J(165)
      COMMON/OVRLAP/MU,NU,MUP,NUP,NONORT,NOVLPS,IROWMU,IROWNU,ICOLMU,
     : ICOLNU,NORTH,IORDER,NCALLS,LMU,LNU,LMUP,LNUP,JMU,JNU,JMUP,JNUP,
     :     IORTH(NORD)
      COMMON/STORE/I(165),I1,I2,I3,I4,I5,I6,I7,I8
*     STORES THE COMMON BLOCK MEDEFN , AND IRHO,ISIG,IRHOP,ISIGP
*
      DO 1 K=1,165
      I(K)=J(K)
    1 CONTINUE
      I1=IRHO
      I2=ISIG
      I3=IRHOP
      I4=ISIGP
      IF(NOVLPS.EQ.0) RETURN
      I5=MU
      I6=NU
      I7=MUP
      I8=NUP
      RETURN
      END
*
*     ------------------------------------------------------------------
*       M E R E S T
*     ------------------------------------------------------------------
*
      SUBROUTINE MEREST(IRHO,ISIG,IRHOP,ISIGP)
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      PARAMETER (NWD=30,NCD=100,NORD=NWD*(NWD-1)/2)
*
      COMMON/MEDEFN/J(165)
      COMMON/OVRLAP/MU,NU,MUP,NUP,NONORT,NOVLPS,IROWMU,IROWNU,ICOLMU,
     : ICOLNU,NORTH,IORDER,NCALLS,LMU,LNU,LMUP,LNUP,JMU,JNU,JMUP,JNUP,
     :     IORTH(NORD)
      COMMON/STORE/I(165),I1,I2,I3,I4,I5,I6,I7,I8
*
*     RESTORES THE COMMON BLOCK MEDEFN, AND IRHO,ISIG,IRHOP,ISIGP
*
      DO 1 K=1,165
      J(K)=I(K)
    1 CONTINUE
      IRHO=I1
      ISIG=I2
      IRHOP=I3
      ISIGP=I4
      IF(NOVLPS.EQ.0) RETURN
      MU=I5
      NU=I6
      MUP=I7
      NUP=I8
      RETURN
      END
*
*     ------------------------------------------------------------------
*       M O D J 2 3
*     ------------------------------------------------------------------
*
      SUBROUTINE MODJ23(KK)
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      PARAMETER (NWD=30,NCD=100,NORD=NWD*(NWD-1)/2)
*
      COMMON/INFORM/IREAD,IWRITE,IOUT,ISC0,ISC1,ISC2,ISC3,JSC(3),IALL
      COMMON/DEBUG/IBUG1,IBUG2,IBUG3,NBUG6,NBUG7,IFULL
      COMMON/MVALUE/M1,M2,M3,M4,M5,M6,M7,M8,M9,M10,M11,M12,M13,M14,M15,
     :     M16,M17,M18,M19,M20
      PARAMETER(KFL1=60,KFL2=12)
      LOGICAL FREE
      COMMON/COUPLE/NJ1S,NJ23S,J1(KFL1),J2(KFL2,3),J3(KFL2,3),FREE(KFL1)
      COMMON/HOLD/J2SPIN(3*KFL2-3),J3SPIN(3*KFL2-3),
     :            J2ANG(3*KFL2),J3ANG(3*KFL2)
      COMMON/OVRLAP/MU,NU,MUP,NUP,NONORT,NOVLPS,IROWMU,IROWNU,ICOLMU,
     : ICOLNU,NORTH,IORDER,NCALLS,LMU,LNU,LMUP,LNUP,JMU,JNU,JMUP,JNUP,
     :     IORTH(NORD)
*
* === MODIFIES THE DIRECT J2 AND J3 ARRAYS FOR EXCHANGE CALL OF NJSYM
*
    7 FORMAT(3H J2,18X,2HJ3)
    8 FORMAT(3I5,I10,2I5)
      GO TO (1,2,1,2,1,2),KK
*
* --- KK=1 OR KK=3 -  SPIN INTEGRALS
*
    1 MK=NJ23S-1
      I1=0
      DO 11 J=1,MK
      DO 12 K=1,3
      I1=I1+1
      J2(J,K)=J2SPIN(I1)
      J3(J,K)=J3SPIN(I1)
   12 CONTINUE
   11 CONTINUE
      IF(KK.EQ.3) GO TO 15
      I1 = 1
      IF (NOVLPS .LT. 2) GO TO 13
      IF (MUP .EQ. NUP) I1=2
   13 IP1 = I1 + 1
      J3(I1,2) = M14
      J3(IP1,2) = M13
      IF (KK.EQ.1) GO TO 3
   15 JSTO=J3(IROWMU,ICOLMU)
      J3(IROWMU,ICOLMU)=J3(IROWNU,ICOLNU)
      J3(IROWNU,ICOLNU)=JSTO
      GO TO 3
*
* --- KK=2 OR KK=4 -  ANGULAR INTEGRALS
*
    2 MK=NJ23S-1
      I1=0
      DO 21 J=1,MK
      DO 22 K=1,3
      I1=I1+1
      J2(J,K)=J2ANG(I1)
      J3(J,K)=J3ANG(I1)
   22 CONTINUE
   21 CONTINUE
      IF(KK.EQ.4) GO TO 18
      I1 = 1
      I2 = 1
      IF(NOVLPS .LT.2) GO TO 14
      IF (MU.EQ.NU .AND. MUP.NE.NUP) I1 = 2
      IF (MU.NE.NU .AND. MUP.EQ.NUP) I2 = 2
   14 J2(I1,1) = M15
      J3(I2,3) = M16
      IF (KK.EQ.2) GO TO 3
   18 IROMU1=IROWMU+1
      IRONU1=IROWNU+1
      JSTO=J3(IROMU1,ICOLMU)
      J3(IROMU1,ICOLMU)=J3(IRONU1,ICOLNU)
      J3(IRONU1,ICOLNU)=JSTO
    3 IF(IBUG2-1) 4,9,9
*
*     PRINT-OUT OF VALUES IN NJSYM  IF IBUG3=1
*
    9 IF(IBUG3-1 ) 5,4,5
    5 WRITE(IWRITE,7)
      DO 6 J=1,MK
      WRITE(IWRITE,8) (J2(J,K),K=1,3),(J3(J,K),K=1,3)
    6 CONTINUE
    4 RETURN
      END
*
*     ------------------------------------------------------------------
*       N O R T B P
*     ------------------------------------------------------------------
*
      SUBROUTINE NORTBP(JA,JB)
      PARAMETER (NWD=30,NCD=100,NORD=NWD*(NWD-1)/2)
*
      DIMENSION ILNO(2),IRNO(2)
      COMMON/INFORM/IREAD,IWRITE,IOUT,ISC0,ISC1,ISC2,ISC3,JSC(3),IALL
      COMMON/DEBUG/IBUG1,IBUG2,IBUG3,NBUG6,NBUG7,IFULL
      COMMON/OVRLAP/MU,NU,MUP,NUP,NONORT,NOVLPS,IROWMU,IROWNU,ICOLMU,
     : ICOLNU,NORTH,IORDER,NCALLS,LMU,LNU,LMUP,LNUP,JMU,JNU,JMUP,JNUP,
     :     IORTH(NORD)
      COMMON/STATES/NCFG,MAXORB,IAJCMP(NWD),LJCOMP(NWD),
     :NJCOMP(NWD),NOCCSH(NCD),NELCSH(5,NCD),NOCORB(5,NCD),J1QNRD(9,NCD)
*
* === Determines, in each of the two configurations JA and JB in a 
*     matrix element, which subshells are non-orthogonal.
*
  101 FORMAT(/63H INCORRECT NON-ORTHOGONALITY SET UP IN THE MATRIX ELEME
     :NT  -  (,I3,3H/V/,I3,1H))
*
      N1=NOCCSH(JA)
      N2=NOCCSH(JB)
      JMU=0
      JNU=0
      JMUP=0
      JNUP=0
      IL=0
      IR=0
*
* --- BEGIN SEARCH FOR NON-ORTHOGONAL SUBSHELLS IN THIS MATRIX ELEMENT
*
      DO 1 I=1,N1
      NI=NOCORB(I,JA)
      DO 2 J=1,N2
      NJ=NOCORB(J,JB)
      IF(NI.EQ.NJ) GO TO 2
      NA=MIN0(NI,NJ)
      NB=MAX0(NI,NJ)
      NC=(NB-2)*(NB-1)/2+NA
      IF(IORTH(NC).NE.1) GO TO 2
      IF (IL .EQ. 0) GO TO 4
      IF (ILNO(IL) .EQ. I) GO TO 14
      IF (IL .EQ. 2) GO TO 100
    4 IL = IL+1
      ILNO(IL) = I
   14 IF ( IR .EQ. 0) GO TO 7
      DO 15 K = 1,IR
      IF (IRNO(K) .EQ. J) GO TO 2
   15 CONTINUE
      IF (IR .EQ. 2) GO TO 100
    7 IR = IR+1
      IRNO(IR) = J
    2 CONTINUE
    1 CONTINUE
      IF(IL.EQ.0) GO TO 8
      IF(IR.EQ.1) GO TO 11
      IF(IRNO(1).LE.IRNO(2)) GO TO 11
      ISTO=IRNO(1)
      IRNO(1)=IRNO(2)
      IRNO(2)=ISTO
   11 JMU=ILNO(1)
    3 IF(IL.EQ.1) GO TO 5
      JNU=ILNO(2)
    5 JMUP=IRNO(1)
    6 IF (IR .EQ. 1) GO TO 10
      JNUP=IRNO(2)
      GO TO 10
    8 NOVLPS=0
      RETURN
  100 WRITE(IWRITE,101) JA,JB
      STOP
   10 NMU=NOCORB(JMU,JA)
      NMUP=NOCORB(JMUP,JB)
      LMU=LJCOMP(NMU)
      LMUP=LJCOMP(NMUP)
      IF (JNU .EQ. 0 .AND. JNUP .EQ. 0) GO TO 9
      IF (JNU .EQ. 0) JNU = JMU
      IF (JNUP .EQ. 0) JNUP = JMUP
      NNU=NOCORB(JNU,JA)
      NNUP=NOCORB(JNUP,JB)
      LNU=LJCOMP(NNU)
      LNUP=LJCOMP(NNUP)
      NOVLPS=2
      RETURN
    9 NOVLPS=1
      RETURN
      END
*
*     ------------------------------------------------------------------
*       N S C R A P
*     ------------------------------------------------------------------
*
      SUBROUTINE NSCRAP(IX,IRS,IR1)
      PARAMETER(KFL2=12)
      DIMENSION IX(KFL2,3)
*
* === Create a blank row in array IX after row (IR1 -1)
*
      IR2=IRS+1-IR1
      DO 1 I=1,IR2
      II=IRS+1-I
      I1=II+1
      DO 2 K=1,3
      IX(I1,K)=IX(II,K)
    2 CONTINUE
    1 CONTINUE
      IRS=IRS+1
      RETURN
      END
*
*     ------------------------------------------------------------------
*       O R T H O G
*     ------------------------------------------------------------------
*
      SUBROUTINE ORTHOG(LET)
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      COMMON/INFORM/IREAD,IWRITE,IOUT,ISC0,ISC1,ISC2,ISC3,JSC(3),IALL
      COMMON/DEBUG/IBUG1,IBUG2,IBUG3,NBUG6,NBUG7,IFULL
      COMMON/MEDEFN/IHSH,NJ(10),LJ(10),NOSH1(10),NOSH2(10),J1QN1(19,3),
     :     J1QN2(19,3),IJFUL(10)
      COMMON/SIGNF /SIGNFA
*
*     THIS SUBROUTINE CHECKS FOR POSSIBLE ORTHOGONALITY DUE TO
*     COUPLING DIFFERENCES OR UNEVEN PARITY
*
  101 FORMAT(37H DIFFERING RESULTANT ANGULAR MOMENTUM)
  102 FORMAT(52H ORTHOGONALITY IN COUPLING SCHEMES OF CONFIGURATIONS)
  103 FORMAT(59H THE TWO CONFIGURATIONS HAVE DIFFERING NUMBERS OF ELECTR
     :ONS)
  104 FORMAT(51H THE TWO CONFIGURATIONS HAVE DIFFERING TOTAL PARITY)
*
* --- DO PSI AND PSIP CONTAIN THE SAME NUMBERS OF ELECTRONS
*     DO PSI AND PSIP HAVE THE SAME TOTAL PARITY
*
      N5=0
      N6=0
      N7=0
      DO 20 I=1,IHSH
      L1=LJ(I)
      L2=NOSH1(I)
      L3=NOSH2(I)
      N5=N5+L2
      N6=N6+L3
      N7=N7+L1*(L2-L3)
   20 CONTINUE
*
*     CHECK ON NUMBER OF ELECTRONS
*
      IF (N5-N6) 21,22,21
   21 CONTINUE
   28 WRITE(IWRITE,103)
      GO TO 11
*
*     CHECK ON PARITY
*
   22 IF(N7-N7/2*2) 23,24,23
   23 CONTINUE
   25 WRITE(IWRITE,104)
      GO TO 11
   24 N72=N7/2
      SIGNFA=1.D0
      IF( (N72-(N72/2)*2).NE.0 ) SIGNFA=-1.D0
      N1=2*IHSH-1
      N2=IHSH+1
      N3=IHSH-1
      N4=IHSH-2
*
* --- IS THE FINAL STATE THE SAME FOR PSI AND PSIP
*
      DO 1 K=2,3
      IF(J1QN1(N1,K)-J1QN2(N1,K))2,1,2
    1 CONTINUE
      GO TO 3
    2 CONTINUE
*  26 WRITE(IWRITE,101)
*
* --- THE TWO CONFIGURATIONS WILL HAVE ZERO HAMILTONIAN MATRIX ELEMENT
*
   11 LET=0
      RETURN
*
* --- COUPLING ORTHOGONALITY TEST FOR FIRST TWO SHELLS
*
    3 CONTINUE
   12 LET=1
      RETURN
      END

*     ------------------------------------------------------------------
*       O U T L S
*     ------------------------------------------------------------------
*
      SUBROUTINE OUTLS
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      PARAMETER (NWD=30,NCD=100)
      COMMON /CLOSED/B1ELC(4),NCLOSD,IAJCLD(NWD),LJCLSD(NWD),IBK
      COMMON/DEBUG/IBUG1,IBUG2,IBUG3,NBUG6,NBUG7,IFULL
      COMMON /FOUT/NOV(2),IOVLAP(10,2),NF,NG,NR,NL,NZ,NN,NV,NS,IFLAG,NIJ
      COMMON/INFORM/IREAD,IWRITE,IOUT,ISC(4),JSC(3),IALL
      COMMON/STATES/NCFG,MAXORB,IAJCMP(NWD),LJCOMP(NWD),
     :NJCOMP(NWD),NOCCSH(NCD),NELCSH(5,NCD),NOCORB(5,NCD),J1QNRD(9,NCD)
      PARAMETER (NCDIM=3000,LSTACK=30)
      DIMENSION C(NCDIM),IPACK(NCDIM),JA(NCDIM),JB(NCDIM),IPTR(NCDIM),
     :          IPT(NCDIM)
      INTEGER NCOUNT(4),ISTACK(LSTACK),II(4),IEL(4)
      EQUIVALENCE (NCOUNT(1),NF),(II(1),I1),(II(2),I2),(II(3),I3),
     :            (II(4),I4)
      CHARACTER*1 INT(8)
      DATA INT/'F','G','R','L','Z','N','V','S'/
*
*     Formats for Integrals
*
   10 FORMAT(1X,A1,I2,'(',A3,',',A3,')',I5)
   11 FORMAT(1X,A1,I2,'(',2A3,',',2A3,')',I5)
   12 FORMAT(1X,A1,2X,'(',A3,',',A3,')',I5)
   13 FORMAT(1X,2('O',I2,'<',A3,'|',A3,'>':))
*     
*     Format for Coefficients
*
   20 FORMAT(F14.8,A1,3I3)
*
*     Format for Coefficient terminator
*
   30 FORMAT(14X,'*',16X,' ')
   31 FORMAT(14X,'*',22X,' ')
   32 FORMAT(14X,'*',14X,' ')
   34 FORMAT(1X,'*')
  103 FORMAT(4(1X,'*'/14X,'*',14X,' '/))
*
* --- TERMINATE  INTEGERAL LISTS and Rewind
*
      DO 1 J = 1,4
         ENDFILE(UNIT=ISC(J))
         REWIND(UNIT=ISC(J))
    1 CONTINUE
*
*     Test if current dimensions are big enough
*
      N = MAX(NF,NG,NR,NL,NZ,NN,NV,NS)
      IF (N .GT. NCDIM) THEN
         WRITE(0,'(A/A,I10)') ' NCDIM dimension in OUT too small',
     :    ' Must be increased to at least',N
         STOP 1
      END IF
*
*===  Begin processing data
*
      NINT = 0
      DO 100 ICASE = 1,4
         N = NCOUNT(ICASE)
         IF (ICASE .LE. 2) THEN
            DO 102 J = 1,N
              READ(ISC(ICASE)) C(J),IPACK(J),JA(J),JB(J)
  102       CONTINUE
         ELSE
            DO 104 J = 1,N
              READ(ISC(ICASE)) C(J),IPACK(J),JA(J),JB(J),IPTR(J)
  104       CONTINUE
         END IF
         CALL QSORT(N,IPACK,IPT,ISTACK,LSTACK,IERR)
         IF (IERR .EQ. 1) THEN
            WRITE(0,*) ' Stack dimension not large enough for sort',
     :           'CASE = ',icase, 'N = ',n
            STOP 1
         END IF
*
*        Output the list of integrals with pointers to the data
*
         LAST = 0
  110    J = LAST +1
         LAST = J
         IF (J .LE. N) THEN
*
*          Unpack electron data
*
           K = IPACK(IPT(J))
           I4 = MOD(K,64)
           K = K/64
           IF (ICASE.LE.2 .OR. ICASE.EQ.4 .OR. ICASE.EQ.5) THEN
             I2 = MOD(K,64)
             K = K/64
           ELSE
              I3 = MOD(K,64)
              K = K/64
              I2 = MOD(K,64)
              K = K/64
              I1 = MOD(K,64)
              K = K/64
           END IF
*
*           Find  last item in the list with this integral
*           
  120      LAST = LAST + 1
           IF (LAST .LE. N) THEN
             IF (IPACK(IPT(J)) .EQ. IPACK(IPT(LAST))) GO TO 120
           END IF
           LAST = LAST -1
           NINT = NINT + 1
           IF (ICASE .LE. 2) THEN
             WRITE(IOUT,10) INT(ICASE),K,IAJCMP(I2),IAJCMP(I4),LAST
           ELSE IF (ICASE .EQ. 4) THEN
             WRITE(IOUT,12) INT(ICASE),IAJCMP(I2),IAJCMP(I4),LAST
           ELSE
             DO 140 J = 1,4
               IF (II(J) .LT. 32) THEN
                 IEL(J) = IAJCMP(II(J))
               ELSE
                 IEL(J) = IAJCLD(64-II(J))
               END IF
  140        CONTINUE
             WRITE(IOUT,11) INT(ICASE),K,(IEL(J),J=1,4),LAST
           END IF
           GO TO 110
        END IF
        WRITE(IOUT,34) 
*
*             Write out the data for the integrals
*
        DO 150 J = 1,N
          K = IPT(J)
          IF (ICASE .LE. 2) THEN
            WRITE(IOUT,20) C(K),INT(ICASE),JA(K),JB(K)
          ELSE IF (ICASE .LE. 4) THEN
            WRITE(IOUT,20) C(K),INT(ICASE),JA(K),JB(K),IPTR(K)
          END IF
  150   CONTINUE
        IF (ICASE .LE. 2) THEN
           WRITE(IOUT,30) 
        ELSE IF (ICASE .EQ. 3) THEN
           WRITE(IOUT,31)
        ELSE
           WRITE(IOUT,32)
        END IF
        IF (ICASE .EQ. 2) THEN
           M = NOV(1)
           DO 160 J = 1,M
              K = IOVLAP(J,1)
              I2 = MOD(K,64)
              K = k/64
              I1 = MOD(K,64)
              K = K/64
              WRITE(IOUT,13) K,IAJCMP(I1),IAJCMP(I2)
  160      CONTINUE
           WRITE(IOUT,34) 
           M = NOV(2)
           DO 170 J = 1,M
             K = IOVLAP(J,2)
             I4 = MOD(K,64)
             K = K/64
             I3 = MOD(K,64)
             K = K/64
             I2 = MOD(K,64)
             K = K/64
             I1 = MOD(K,64)
             K = K/64
             K1 = MOD(K,16)
             K = K/16
             WRITE(IOUT,13) K,IAJCMP(I1),IAJCMP(I2),
     :                            K1,IAJCMP(I3),IAJCMP(I4)
  170      CONTINUE
           WRITE(IOUT,34)
        END IF
        CLOSE(UNIT=ISC(ICASE))
  100 CONTINUE
      WRITE(IOUT,103)
      WRITE(0,*) 'The total number of integrals =',NINT+NOV(1)+NOV(2)
      RETURN
      END
*-----------------------------------------------------------------------
*        O V L S E T
*-----------------------------------------------------------------------
*
      SUBROUTINE OVLSET(N,I1,I2,IOV,I3,I4,JOV,IPTR)
   
      COMMON /FOUT/NOV(2),IOVLAP(10,2),NF,NG,NR,NL,NZ,NN,NV,NS,IFLAG,NIJ

      IF ( N .EQ. 0) THEN
         IPTR = 0
         RETURN
      ELSE
         IF (I1 .GT. I2) THEN
            IT = I1
            I1 = I2
            I2 = IT
         END IF
         IF (N .EQ. 1) THEN
            IPACK = (IOV*64 + I1)*64 +I2
         ELSE
            IF (I3 .GT. I4) THEN
               IT = I3
               I3 = I4
               I4 = IT
            END IF
            IF (I1 .GT. I3) THEN
               IT = I1
               I1 = I3
               I3 = I1
               IT = I2
               I2 = I4
               I4 = I2
               IT = IOV
               IOV = JOV
               JOV = IT
            END IF
            IPACK = ((((IOV*16+JOV)*64+I1)*64+I2)*64+I3)*64+I4
         END IF
      END IF

      M = NOV(N)
      DO 10 I = 1,M
         IF ( IOVLAP(I,N) .EQ. IPACK) THEN
            IF (N .EQ. 1) THEN
               IPTR = I
            ELSE
               IPTR = -I
            END IF
            RETURN
         END IF
10    CONTINUE
      M = M + 1
      NOV(N) = M
      IOVLAP(M,N) = IPACK
      IF (N .EQ. 1) THEN
         IPTR = M
      ELSE
         IPTR = -M
      END IF
      END
*
*     ------------------------------------------------------------------
*       P I K U P 1
*     ------------------------------------------------------------------
*
      SUBROUTINE PIKUP1(IX,IRS,MU,NU,IR1,IR2,IC1,IC2)
      PARAMETER(KFL2=12)
      DIMENSION IX(KFL2,3)
*
* === Locates the position in array IX of MU and NU
*
      IR1=0
      IR2=0
      DO 1 I=1,IRS
      DO 2 K=1,2
      IA=IX(I,K)
      IF(IA.EQ.MU) GO TO 3
      IF(IA.EQ.NU) GO TO 4
      GO TO 2
    3 IR1=I
      IC1=K
      IF(MU.EQ.NU) RETURN
      GO TO 2
    4 IR2=I
      IC2=K
    2 CONTINUE
    1 CONTINUE
      RETURN
      END
*
*     ------------------------------------------------------------------
*       P I K U P 2
*     ------------------------------------------------------------------
*
      SUBROUTINE PIKUP2(IX,IRS,IR1,IMU,IR2,IC2)
      PARAMETER(KFL2=12)
      DIMENSION IX(KFL2,3)
*
* === Locates IMU in IX
*
      IR2=0
      IR3=IR1+1
      DO 1 I=IR3,IRS
      DO 2 K=1,2
      IA=IX(I,K)
      IF(IA.NE.IMU) GO TO 2
      IR2=I
      IC2=K
      RETURN
    2 CONTINUE
    1 CONTINUE
      RETURN
      END
*
*     ------------------------------------------------------------------
*       P R N T W T
*     ------------------------------------------------------------------
*
      SUBROUTINE PRNTWT(IRHO,ISIG,IRHOP,ISIGP)
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      PARAMETER (NWD=30,NCD=100,NORD=NWD*(NWD-1)/2)
*
      COMMON/OVRINT/IOVEL(2),IOVER(2),IOVEP(2)
      COMMON/INFORM/IREAD,IWRITE,IOUT,ISC0,ISC1,ISC2,ISC3,JSC(3),IALL
      COMMON/DEBUG/IBUG1,IBUG2,IBUG3,NBUG6,NBUG7,IFULL
      COMMON/STATES/NCFG,MAXORB,IAJCMP(NWD),LJCOMP(NWD),
     :NJCOMP(NWD),NOCCSH(NCD),NELCSH(5,NCD),NOCORB(5,NCD),J1QNRD(9,NCD)
      COMMON/ENAV/NINTS,KVALUE(15),COEFCT(15)
      COMMON/MEDEFN/IHSH,NJ(10),LJ(10),NOSH1(10),NOSH2(10),J1QN(19,3,2)
     :,IJFUL(10)
      COMMON/DIAGNL/IDIAG,JA,JB
      COMMON/MVALUE/M1,M2,M3,M4,M5,M6,M7,M8,M9,M10,M11,M12,M13,M14,M15,
     :     M16,M17,M18,M19,M20
      COMMON/XATION/KD1,KD2,KE1,KE2,MULTD,MULTE,
     :     AMULT1(19),AMULT2(19),BMULT1(19),BMULT2(19)
      COMMON/NJLJ/NRHO,LRHO,NSIG,LSIG,NRHOP,LRHOP,NSIGP,LSIGP
      COMMON/SIGNF /SIGNFA
      COMMON/OVRLAP/MU,NU,MUP,NUP,NONORT,NOVLPS,IROWMU,IROWNU,ICOLMU,
     : ICOLNU,NORTH,IORDER,NCALLS,LMU,LNU,LMUP,LNUP,JMU,JNU,JMUP,JNUP,
     :     IORTH(NORD)
*
* --- PRINTS OUT THE COEFFICIENTS OF SLATER INTEGRALS
*
    1 FORMAT(//23H INTERACTING SHELLS ARE,6X,6H RHO =,A3,6X,6H SIG =,A3,
     :6X,7H RHOP =,A3,6X,7H SIGP =,A3//)
    2 FORMAT(3F14.6,4X,1HF,I2,1H(,A3,1H,,A3,1H))
    3 FORMAT(3F14.6,4X,1HG,I2,1H(,A3,1H,,A3,1H))
    4 FORMAT(3F14.6,4X,1HR,I2,1H(,A3,1H,,A3,1H;,A3,1H,,A3,1H))
  104 FORMAT(3F14.6,4X,1HR,I2,1H(,A3,1H,,A3,1H;,A3,1H,,A3,1H),
     :1X,1H<,A3,1H|,A3,3H>**,I2)
  204 FORMAT(3F14.6,4X,1HR,I2,1H(,A3,1H,,A3,1H;,A3,1H,,A3,1H),
     :1X,1H<,A3,1H|,A3,3H>**,I2,2H <,A3,1H|,A3,3H>**,I2)
  100 FORMAT(F14.8,1HF,I2,1H(,A3,I3,1H,,A3,I3,1H))
  101 FORMAT(F14.8,1HG,I2,1H(,A3,I3,1H,,A3,I3,1H))
  102 FORMAT(F14.8,1HR,I2,1H(,2A3,I3,1H,,2A3,I3,1H):
     :       2(1H<,A3,1H|,A3,1H>,I2:))
      JRHO=IJFUL(IRHO)
      JSIG=IJFUL(ISIG)
      JRHOP=IJFUL(IRHOP)
      JSIGP=IJFUL(ISIGP)
   42 IF(NOVLPS.EQ.0) GO TO 40
      JMU = IJFUL(MU)
      JMUP = IJFUL(MUP)
   43 IF (NOVLPS .EQ. 1) GO TO 40
      JNU = IJFUL(NU)
      JNUP = IJFUL(NUP)
*
* --- DETERMINE THE AVERAGE ENERGY CONTRIBUTIONS IF  IDIAG  IS NON-ZERO
*
   40 IF(KD1.GT.KD2.AND.KE1.GT.KE2) GO TO 41
      IF(IDIAG.EQ.0) GO TO 50
      LA=LJ(IRHO)
      LB=LJ(ISIG)
      IF(M1.EQ.0) GO TO 51
      IEQUIV=2
      NMULT=NOSH1(IRHO)*NOSH1(ISIG)
      GO TO 52
   51 IEQUIV=1
      NA=NOSH1(IRHO)
      NMULT=NA*(NA-1)/2
*
*     CALCULATE THE INTERACTION ENERGY
*
   52 CALL INTACT(LA,LB,IEQUIV)
      INTS=1
   50 IF(IBUG2-1)  6,7,7
    7 WRITE(IWRITE,1) IAJCMP(JRHO),IAJCMP(JSIG),IAJCMP(JRHOP),
     :                IAJCMP(JSIGP)
    6 CONTINUE
* --- DIRECT INTEGRALS
*
 8    IF(KD1.GT.KD2) GO TO 20
      DO 11 JK1=KD1,KD2,2
      K=JK1-1
      A=AMULT1(JK1)*SIGNFA
      IF(NONORT.NE.0) A2=AMULT2(JK1)*SIGNFA
*
* --- DIVIDE THE WEIGHTS INTO AVERAGE ENERGY AND NON-AVERAGE ENERGY
*     PARTS
*
      IF(IDIAG.NE.0) GO TO 53
*
*     NON-DIAGONAL MATRIX ELEMENT
*
      X=0.D0
      Y=A
      IF(NONORT.NE.0) Y2=A2
      IF(NOVLPS.NE.0) GO TO 13
      IF(M19.EQ.0.AND.M20.EQ.0) GO TO 15
*     IF((M1+M2).EQ.0) GO TO 16
      GO TO 13
*
*     DIAGONAL MATRIX ELEMENT.  F0 TERM IS THE ONLY ONE WITH K=0
*
   53 IF(K.NE.0) GO TO 57
      X=DFLOAT(NMULT)
      Y=A-X
      GO TO 15
*
*     OTHER  FK  INTEGRALS,  ONLY OCCUR IF  RHO=SIG
*
   57 IF(IEQUIV.EQ.1.AND.INTS.LE.NINTS) GO TO 58
   59 X=0.D0
      Y=A
      GO TO 15
   58 IF(K.NE.KVALUE(INTS)) GO TO 59
      X=NMULT*COEFCT(INTS)
      Y=A-X
      INTS=INTS+1
   15 IF((DABS(A).GE.1.D-10.OR.DABS(X).GE.1.D-10) .AND. IFULL .NE. 0)
     :   WRITE(IWRITE,2)A,X,Y,K,IAJCMP(JRHO),IAJCMP(JSIG)
      IF (IALL .NE. 0) Y = A
      IF (IALL.EQ.0 .AND. JA.NE.JB) Y=Y*2.D0
      IF(DABS(Y).GT.1.D-10) 
     : CALL SAVE(1,Y,K,0,JRHO,0,JSIG,JA,JB,0)
      GO TO 11
   16 IF(DABS(A).GE.1.D-10 .AND. IFULL .NE. 0)
     :   WRITE(IWRITE,3) A,X,Y,K,IAJCMP(JRHO),IAJCMP(JRHOP)
      IF (IALL .NE. 0) Y = A
      IF (IALL.EQ.0 .AND. JA.NE.JB) Y=Y*2.D0
      IF(DABS(Y).GT.1.D-10) 
     :  CALL SAVE(2,Y,K,0,JRHO,0,JRHOP,JA,JB,0)
      GO TO 11
   13 NOVLP1=NOVLPS+1
      GO TO (213,313,413),NOVLP1
  213 IF(DABS(A).GE.1.D-10 .AND.IFULL.NE.0) WRITE(IWRITE,4) A,X,Y,K,
     : IAJCMP(JRHO),IAJCMP(JSIG),IAJCMP(JRHOP),IAJCMP(JSIGP)
      IF (IALL .NE. 0) Y = A
      IF (IALL.EQ.0 .AND. JA.NE.JB) Y = Y+Y
      IF(DABS(Y).GT.1.D-10) 
     :  CALL SAVE(3,Y,K,JRHO,JSIG,JRHOP,JSIGP,JA,JB,0)
      GO TO 11
  313 IF(DABS(A).GE.1.D-10 .AND. IFULL .NE. 0)
     :  WRITE(IWRITE,104) A,X,Y,K,IAJCMP(JRHO),IAJCMP(JSIG),
     : IAJCMP(JRHOP),IAJCMP(JSIGP),IAJCMP(JMU),IAJCMP(JMUP),IOVEP(1)
      IF (IALL .NE. 0) Y = A
      IF (IALL.EQ.0 .AND. JA.NE.JB) Y = Y+Y
      IF(DABS(Y).GT.1.D-10) THEN
        CALL OVLSET(1,JMU,JMUP,IOVEP(1),0,0,0,IPTR)
        CALL SAVE(3,Y,K,JRHO,JSIG,JRHOP,JSIGP,JA,JB,IPTR)
      END IF
      GO TO 11
  413 IF(DABS(A).GE.1.D-10 .AND.IFULL.NE.0) WRITE(IWRITE,204) A,X,Y,K,
     : IAJCMP(JRHO),IAJCMP(JSIG),IAJCMP(JRHOP),IAJCMP(JSIGP),IAJCMP(JMU)
     : ,IAJCMP(JMUP),IOVEP(1),IAJCMP(JNU),IAJCMP(JNUP),IOVEP(2)
      IF (IALL .NE. 0) Y = A
      IF (IALL.EQ.0 .AND. JA.NE.JB) Y = Y+Y
      IF(DABS(Y).GT.1.D-10) THEN
        CALL OVLSET(2,JMU,JMUP,IOVEP(1),JNU,JNUP,IOVEP(2),IPTR)
        CALL SAVE(3,Y,K,JRHO,JSIG,JRHOP,JSIGP,JA,JB,IPTR)
      END IF
      IF(NONORT.EQ.0 ) GO TO 11
      IF (DABS(A2).GE.1.D-10.AND.IFULL.NE.0) WRITE(IWRITE,204)A2,X,Y2,K,
     : IAJCMP(JRHO),IAJCMP(JSIG),IAJCMP(JRHOP),IAJCMP(JSIGP),IAJCMP(JNU)
     : ,IAJCMP(JMUP),IOVEP(2),IAJCMP(JMU),IAJCMP(JNUP),IOVEP(1)
      IF ( IALL .NE. 0) Y2 = A2
      IF (IALL.EQ.0 .AND. JA.NE.JB) Y2 = Y2+Y2
      IF(DABS(Y2).GT.1.D-10) THEN
        CALL OVLSET(2,JNU,JMUP,IOVEP(2),JMU,JNUP,IOVEP(1),IPTR)
        CALL SAVE(3,Y2,K,JRHO,JSIG,JRHOP,JSIGP,JA,JB,IPTR)
      END IF
   11 CONTINUE
*
* --- EXCHANGE INTEGRALS
*
 20   IF(KE1.GT.KE2) GO TO 41
      DO 21 JK1=KE1,KE2,2
      K=JK1-1
      B=BMULT1(JK1)*SIGNFA
      IF(NONORT.NE.0) B2=BMULT2(JK1)*SIGNFA
*
* --- DIVIDE THE WEIGHTS INTO AVERAGE ENERGY AND NON-AVERAGE ENERGY
*     PARTS
*
      IF(IDIAG.NE.0) GO TO 60
*
*     NON-DIAGONAL MATRIX ELEMENT
*
      X=0.D0
      Y=B
      IF (NONORT .NE. 0) Y2 = B2
      IF(M19.EQ.0.AND.M20.EQ.0) GO TO 25
      IF(M1 .EQ.0.AND.M2 .EQ.0) GO TO 31
      GO TO 23
*
*     DIAGONAL MATRIX ELEMENT
*
   60 IF(INTS.LE.NINTS) GO TO 61
   62 X=0.D0
      Y=B
      GO TO 25
   61 IF(K.NE.KVALUE(INTS)) GO TO 62
      X=NMULT*COEFCT(INTS)
      Y=B-X
      INTS=INTS+1
   25 IF(DABS(B)+DABS(X) .GE. 1.D-10 .AND. IFULL .NE. 0)
     : WRITE(IWRITE,3) B,X,Y,K,IAJCMP(JRHO),IAJCMP(JSIG)
      IF (IALL .NE. 0) Y = B
      IF (IALL.EQ.0 .AND. JA.NE.JB) Y = Y+Y
      IF (DABS(Y) .GE. 1.D-10) 
     :  CALL SAVE(2,Y,K,0,JRHO,0,JSIG,JA,JB,0)
      GO TO 21
   31 IF(DABS(B)+DABS(X) .GE. 1.D-10 .AND. IFULL .NE. 0)
     :  WRITE(IWRITE,3) B,X,Y,K,IAJCMP(JRHO),IAJCMP(JSIGP)
      IF (IALL .NE. 0) Y = B
      IF (IALL.EQ.0 .AND. JA.NE.JB) Y = Y+Y
      IF (DABS(Y) .GE. 1.D-10) 
     :  CALL SAVE(2,Y,K,0,JRHO,0,JSIGP,JA,JB,0)
      GO TO 21
   23 NOVLP1=NOVLPS+1
      GO TO (223,323,423),NOVLP1
  223 IF(DABS(B).GE.1.D-10 .AND. IFULL.NE.0) WRITE(IWRITE,4) B,X,Y,K,
     : IAJCMP(JRHO),IAJCMP(JSIG),IAJCMP(JSIGP),IAJCMP(JRHOP)
      IF (IALL .NE. 0) Y = B
      IF (IALL.EQ.0 .AND. JA.NE.JB) Y = Y+Y
      IF (DABS(Y) .GE. 1.D-10)
     :   CALL SAVE(3,Y,K,JRHO,JSIG,JSIGP,JRHOP,JA,JB,0)
      GO TO 21
  323 IF (DABS(B) .GE. 1.D-10 .AND. IFULL .NE. 0)
     :   WRITE(IWRITE,104) B,X,Y,K,IAJCMP(JRHO),IAJCMP(JSIG),
     :                             IAJCMP(JSIGP),IAJCMP(JRHOP),
     :                             IAJCMP(JMU),IAJCMP(JMUP),IOVEP(1)
      IF (IALL .NE. 0) Y = B
      IF (IALL.EQ.0 .AND. JA.NE.JB) Y = Y+Y
      IF (DABS(Y) .GE. 1.D-10) THEN
         CALL OVLSET(1,JMU,JMUP,IOVEP(1),0,0,0,IPTR)
         CALL SAVE(3,Y,K,JRHO,JSIG,JSIGP,JRHOP,JA,JB,IPTR)
      END IF
      GO TO 21
  423 IF(DABS(B).GE.1.D-10.AND.IFULL.NE.0) WRITE(IWRITE,204) B,X,Y,K,
     : IAJCMP(JRHO),IAJCMP(JSIG),IAJCMP(JSIGP),IAJCMP(JRHOP),IAJCMP(JMU)
     : ,IAJCMP(JMUP),IOVEP(1),IAJCMP(JNU),IAJCMP(JNUP),IOVEP(2)
      IF (IALL .NE. 0) Y = B
      IF (IALL.EQ.0 .AND. JA.NE.JB) Y = Y+Y
      IF (DABS(Y) .GE. 1.D-10) THEN
         CALL OVLSET(2,JMU,JMUP,IOVEP(1),JNU,JNUP,IOVEP(2),IPTR)
         CALL SAVE(3,Y,K,JRHO,JSIG,JSIGP,JRHOP,JA,JB,IPTR)
      END IF
      IF(NONORT.EQ.0) GO TO 21
      IF (DABS(B2).GE.1.D-10.AND.IFULL.NE.0) WRITE(IWRITE,204)B2,X,Y2,K,
     : IAJCMP(JRHO),IAJCMP(JSIG),IAJCMP(JSIGP),IAJCMP(JRHOP),IAJCMP(JNU)
     : ,IAJCMP(JMUP),IOVEP(2),IAJCMP(JMU),IAJCMP(JNUP),IOVEP(1)
      IF ( IALL .NE. 0) Y2 = B2
      IF (IALL.EQ.0 .AND. JA.NE.JB) Y2 = Y2+Y2
      IF (DABS(Y2) .GE. 1.D-10) THEN
         CALL OVLSET(2,JNU,JMUP,IOVEP(2),JMU,JNUP,IOVEP(1),IPTR)
         CALL SAVE(3,Y2,K,JRHO,JSIG,JSIGP,JRHOP,JA,JB,IPTR)
      END IF
   21 CONTINUE
   41 RETURN
      END
*
*     ------------------------------------------------------------------
*       R E D U C E
*     ------------------------------------------------------------------
*
      SUBROUTINE REDUCE(IRHO,ISIG,IRHOP,ISIGP,LESSEN)
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      PARAMETER (NWD=30,NCD=100,NORD=NWD*(NWD-1)/2)
*
      DIMENSION LEAVE(10)
      COMMON/INFORM/IREAD,IWRITE,IOUT,ISC0,ISC1,ISC2,ISC3,JSC(3),IALL
      COMMON/DEBUG/IBUG1,IBUG2,IBUG3,NBUG6,NBUG7,IFULL
      COMMON/MEDEFN/IHSH,NJ(10),LJ(10),NOSH1(10),NOSH2(10),J1QN1(19,3),
     :     J1QN2(19,3),IJFUL(10)
      COMMON/OVRLAP/MU,NU,MUP,NUP,NONORT,NOVLPS,IROWMU,IROWNU,ICOLMU,
     : ICOLNU,NORTH,IORDER,NCALLS,LMU,LNU,LMUP,LNUP,JMU,JNU,JMUP,JNUP,
     :     IORTH(NORD)
*
*     THIS SUBROUTINE REMOVES SPECTATOR SINGLET  S  SHELLS WHICH HAVE
*     NO EFFECT IN ANGULAR OR SPIN INTEGRALS
*
*     LMIN  INITIALLY SET LARGE
*
      LMIN=99
      ICOUNT=0
      DO 1 I=1,IHSH
*
*     NO INTERACTING SHELL MAY BE REMOVED
*
      IF(I.EQ.IRHO.OR.I.EQ.ISIG.OR.I.EQ.IRHOP.OR.I.EQ.ISIGP) GO TO 2
      IF(NOVLPS.EQ.0) GO TO 11
      IF(I.EQ.MU.OR.I.EQ.NU.OR.I.EQ.MUP.OR.I.EQ.NUP) GO TO 2
*
*     IF A SPECTATOR SHELL HAS   SINGLET  S   COUPLING ON BOTH SIDES OF
*     THE MATRIX ELEMENT,  IT MAY, IN GENERAL,  BE REMOVED, AS IT HAS NO
*     EFFECT IN FANO
*
   11 IF(J1QN1(I,1).EQ.0.AND.J1QN2(I,1).EQ.0) GO TO 7
    2 ICOUNT=ICOUNT+1
      LEAVE(ICOUNT)=I
      GO TO 1
    7 IF(LJ(I).GE.LMIN) GO TO 1
      LMIN=LJ(I)
      ILMIN=I
    1 CONTINUE
      IF(ICOUNT.EQ.IHSH) GO TO 8
*
*     IF A CHANGE IN THE COMMON BLOCK  MEDEFN  IS TO BE MADE,
*     ITS PRESENT SITUATION MUST BE PRESERVED BY A CALL OF MEKEEP
*
      CALL MEKEEP(IRHO,ISIG,IRHOP,ISIGP)
*
*     IF ONLY ONE SHELL WOULD BE LEFT IN THIS WAY, THE ONE, DESTINED
*     FOR REMOVAL, WITH THE LOWEST L-VALUE MUST BE RETAINED TO DEFINE A
*     COUPLING
*
      IF(ICOUNT.EQ.1) GO TO 10
*
* --- MODIFY THE COMMON BLOCK  MEDEFN
13    CONTINUE
      DO 3 I=1,ICOUNT
      J=LEAVE(I)
      IF(J.EQ.IRHO) IRHO=I
      IF(J.EQ.ISIG) ISIG=I
      IF(J.EQ.IRHOP) IRHOP=I
      IF(J.EQ.ISIGP) ISIGP=I
      IF(NOVLPS.EQ.0) GO TO 12
      IF(J.EQ.MU) MU=I
      IF(J.EQ.NU) NU=I
      IF(J.EQ.MUP) MUP=I
      IF(J.EQ.NUP) NUP=I
   12 NJ(I)=NJ(J)
      LJ(I)=LJ(J)
      NOSH1(I)=NOSH1(J)
      NOSH2(I)=NOSH2(J)
      DO 4 K=1,3
      J1QN1(I,K)=J1QN1(J,K)
      J1QN2(I,K)=J1QN2(J,K)
    4 CONTINUE
    3 CONTINUE
      ISUBH=IHSH-1
      DO 5 I=2,ICOUNT
      J=LEAVE(I)
      II=ICOUNT+I-1
      IJ=ISUBH+J
      DO 6 K=1,3
      J1QN1(II,K)=J1QN1(IJ,K)
      J1QN2(II,K)=J1QN2(IJ,K)
    6 CONTINUE
    5 CONTINUE
      IHSH=ICOUNT
      GO TO 20
*
*     THIS SITUATION ONLY OCCURS IF IRHO=ISIG=IRHOP=ISIGP
*
   10 J=LEAVE(1)
      II=MIN0(J,ILMIN)
      LEAVE(1)=II
      LEAVE(2)=J+ILMIN-II
      ICOUNT=2
      GO TO 13
   20 IF(IBUG2.EQ.0) GO TO 9
      I2HSH=IHSH+IHSH-1
      WRITE(IWRITE,35)
      WRITE(IWRITE,40) ((J1QN1(J,K),K=1,3),J=1,I2HSH)
      WRITE(IWRITE,36)
      WRITE(IWRITE,40) ((J1QN2(J,K),K=1,3),J=1,I2HSH)
   35 FORMAT(/35H NEW DEFINITION OF COUPLING SCHEMES/38H FOR THIS SET OF
     :  RHO, SIG, RHOP, SIGP//10X,48H L.H.S. OF HAMILTONIAN MATRIX ELEME
     :NT DEFINED BY)
   36 FORMAT(10X,48H R.H.S. OF HAMILTONIAN MATRIX ELEMENT DEFINED BY)
   40 FORMAT(10X,6H J1QN ,9(I5,2I3))
*
*     LESSEN = 0   IF NO CHANGE IN MEDEFN
*            = 1   OTHERWISE
*
    9 LESSEN=1
      RETURN
    8 LESSEN=0
      RETURN
      END
*
*     ------------------------------------------------------------------
*       R K W T S
*     ------------------------------------------------------------------
*
      SUBROUTINE RKWTS
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      PARAMETER (NWD=30,NCD=100)
*
      COMMON/INFORM/IREAD,IWRITE,IOUT,ISC0,ISC1,ISC2,ISC3,JSC(3),IALL
      COMMON/DEBUG/IBUG1,IBUG2,IBUG3,NBUG6,NBUG7,IFULL
      COMMON/DIAGNL/IDIAG,JA,JB
      COMMON /CLOSED/B1ELC(4),NCLOSD,IAJCLD(NWD),LJCLSD(NWD),IBK
      COMMON/MEDEFN/IHSH,NJ(10),LJ(10),NOSH1(10),NOSH2(10),J1QN1(19,3),
     :     J1QN2(19,3),IJFUL(10)
      COMMON/MVALUE/M1,M2,M3,M4,M5,M6,M7,M8,M9,M10,M11,M12,M13,M14,M15,
     :     M16,M17,M18,M19,M20
      COMMON/REMOVE/ICHOP(10)
      COMMON/OVRLAP/MU,NU,MUP,NUP,NONORT,NOVLPS,IROWMU,IROWNU,ICOLMU,
     : ICOLNU,NORTH,IORDER,NCALLS,LMU,LNU,LMUP,LNUP,JMU,JNU,JMUP,JNUP,
     :     IORTH(NWD*(NWD-1)/2)
*
*     THE MATRIX ELEMENT OF THE TWO-ELECTRON POTENTIAL BETWEEN TWO
*     STATES (LABELLED 1 AND 2) MAY BE EXPRESSED AS A SUM OF WEIGHTED
*     RK (SLATER) INTEGRALS.  THIS SUBROUTINE, TOGETHER WITH THOSE
*     CALLED BY IT, DETERMINES THESE WEIGHTS, WHICH ARISE FROM AN
*     INTEGRATION OVER THE ANGULAR AND SPIN CO-ORDINATES
*     FOR DETAILS, SEE   U. FANO, PHYS. REV.,140,A67,(1965)
*
*     THE =INTERACTING= SHELLS ARE DESIGNATED  IRHO,ISIG,IRHOP,ISIGP.
*     THE FIRST TWO REFER TO THE L.H.S. OF     (PSI/V/PSIP)     , WHILE
*     THE SECOND TWO REFER TO THE R.H.S.  FOR DIAGONAL AND CERTAIN OFF-
*     DIAGONAL MATRIX ELEMENTS, THESE MAY NOT BE UNIQUE, AND EACH
*     POSSIBILITY MUST BE CONSIDERED IN TURN
*     THE CONDITION =IRHO .LE. ISIG ,  IRHOP .LE. ISIGP=  IS TO BE
*     SATISFIED
*
   61 FORMAT(//10X,7H IRHO =,I3,4X,7H ISIG =,I3,4X,8H IRHOP =,I3,3X,8H I
     :SIGP =,I3)
*
* === DETERMINE THE INTERACTING SHELLS AS FAR AS POSSIBLE BY
*     CONSIDERING THE DIFFERENCES BETWEEN PSI AND PSIP
*
      IX=0
      IY=0
      IBACK3=0
      NOVSTO = NOVLPS
      IBACK3 = 0
      IBK = 0
      IRHO=0
      ISIG=0
      IRHOP=0
      ISIGP=0
      IZERO = 0
      IA=0
      IB=0
      IC=0
      ID=0
      IONE = 1
      ITWO = 2
      DO 4 J=1,IHSH
      IF(NOVLPS.EQ.0) GO TO 38
      IF(J.EQ.MU.OR.J.EQ.NU.OR.J.EQ.MUP.OR.J.EQ.NUP) GO TO 4
   38 N=NOSH1(J)-NOSH2(J)
      IF(IABS(N)-2) 5,5,1
    5 IF(N) 7,4,6
    6 IF(N-1) 4,8,9
*
* --- TO SATISFY =IRHO.LE.ISIG=  IRHO IS SET FIRST,  ETC.
*
    8 IF(IRHO) 10,10,11
   10 IRHO = J
      GO TO 12
   11 ISIG=J
   12 IX =IX+1
      GO TO 4
    9 IRHO=J
      IX=IX+2
      GO TO 4
    7 IF(N+1) 13,14,4
   14 IF(IRHOP) 15,15,16
   15 IRHOP = J
      GO TO 17
   16 ISIGP=J
   17 IY=IY+1
      GO TO 4
   13 IRHOP=J
      IY=IY+2
    4 CONTINUE
*
*     IX MEASURES THE TOTAL NUMBER OF ELECTRONS IN EITHER CONFIGURATION
*     WHICH DO NOT OCCUR IN THE OTHER.  THEN IF  IX  IS GREATER THAN 4,
*     ORTHOGONALITY OF THE ORBITALS PREVENTS A NON-ZERO MATRIX ELEMENT.
*     IF  IX  IS LESS THAN 4, THEN WE DIVIDE IX BY 2 AND NOW IX MEASURES
*     THE NUMBER OF ELECTRONS WHICH HAVE BEEN CHANGED IN GOING FROM PSI
*     TO PSIP.  IF NOW IX=0, WE HAVE A DIAGONAL MATRIX ELEMENT.  RHO AND
*     SIG MAY TAKE ON ANY VALUES LESS THAN IHSH.  IF IX=1, ONE INTER-
*     ACTING SHELL ON EACH SIDE IS FIXED, WHILE THE OTHER MAY VARY.  IF
*     IX=2, ALL INTERACTING SHELLS ARE DETERMINED
*
      IF(IX.GT.2.OR.IY.GT.2) RETURN
      IF(NOVLPS.EQ.0) GO TO 101
      MUSTO=MU
      NUSTO=NU
      MUPSTO=MUP
      NUPSTO=NUP
      NOVSTO=NOVLPS
      LMUSTO=LMU
      LNUSTO=LNU
      LMUPST=LMUP
      LNUPST=LNUP
  101 IBAL=IABS(IX-IY)
      IF(IX.NE.IY) GO TO 102
      IF(IX-1) 19,20,21
  102 IF(IX.GT.IY) GO TO 121
      IC=0
      ID=0
*
* --- IX.LT.IY
*
      IF(IBAL.EQ.2) THEN
*
*  IX=0, IY=2
*
         IF(IBACK3.EQ.0) THEN
            IRHO=MU
            ISIG=MU
            IA=2
            IB=0
            IBACK3=1
            IF(MU.EQ.NU.OR.NU.EQ.0) IBACK3=0
            GO TO 21
         ELSEIF (IBACK3.EQ.1) THEN
            IRHO=MU
            ISIG=NU
            IA=1
            IB=1
            IBACK3=2
            GO TO 21
         ELSE
            IRHO=NU
            ISIG=NU
            IA=0
            IB=2
            IBACK3=0
            GO TO 21
         ENDIF
      ELSEIF (IX.EQ.1) THEN
*
*   IX=1, IY=2
*
         IF(IBACK3.EQ.0) THEN
            ISIG=MU
            IA=1
            IB=0
            IBACK3=1
            IRSTO=IRHO
            IF(IRHO.GT.ISIG) THEN
               ISTO=IRHO
               IRHO=ISIG
               ISIG=ISTO
            ENDIF
            IF(MU.EQ.NU.OR.NU.EQ.0) IBACK3=0
            GO TO 21
         ELSE
            IRHO=IRSTO
            ISIG=NU
            IA=0
            IB=1
            IBACK3=0
            GO TO 21
         ENDIF
      ELSE
*
*   IX=0, IY=1
*
         IF(IBACK3.EQ.0) THEN
            IRHO=MU
            IA=1
            IB=0
            IBACK3=1
            IF(MU.EQ.NU.OR.NU.EQ.0) IBACK3=0
            GO TO 20
         ELSE
            IRHO=NU
            IA=0
            IB=1
            IBACK3=0
            IBK = 2
            GO TO 20
         ENDIF
      ENDIF
*
*   IX.GT.IY
*
  121 IA=0
      IB=0
      IF(IBAL.EQ.2) THEN
*
*   IX=2, IY=0
*
         IF(IBACK3.EQ.0) THEN
            IRHOP=MUP
            ISIGP=MUP
            IC=2
            ID=0
            IBACK3=1
            IF(MUP.EQ.NUP.OR.NUP.EQ.0) IBACK3=0
            GO TO 21
         ELSEIF (IBACK3.EQ.1) THEN
            IRHOP=MUP
            ISIGP=NUP
            IC=1
            ID=1
            IBACK3=2
            GO TO 21
         ELSE
            IRHOP=NUP
            ISIGP=NUP
            IC=0
            ID=2
            IBACK3=0
            GO TO 21
         ENDIF
      ELSEIF (IY.EQ.1) THEN
*
*   IX=2, IY=1
*
         IF(IBACK3.EQ.0) THEN
            ISIGP=MUP
            IC=1
            ID=0
            IBACK3=1
            IRPSTO=IRHOP
            IF(IRHOP.GT.ISIGP) THEN
               ISTO=IRHOP
               IRHOP=ISIGP
               ISIGP=ISTO
            ENDIF
            IF(MUP.EQ.NUP.OR.NUP.EQ.0) IBACK3=0
            GO TO 21
         ELSE
            IRHOP=IRPSTO
            ISIGP=NUP
            IC=0
            ID=1
            IBACK3=0
            GO TO 21
         ENDIF
      ELSE
*
*   IX=1, IY=0
*
         IF(IBACK3.EQ.0) THEN
            IRHOP=MUP
            IC=1
            ID=0
            IBACK3=1
            IF(MUP.EQ.NUP.OR.NUP.EQ.0) IBACK3=0
            GO TO 20
         ELSE
            IRHOP=NUP
            IC=0
            ID=1
            IBACK3=0
            IBK=2
            GO TO 20
         ENDIF
      ENDIF
*
* === UNIQUE SPECIFICATION OF INTERACTING SHELLS
*
   21 IF(ISIG) 22,23,22
   23 ISIG=IRHO
   22 IF(ISIGP) 24,25,24
   25 ISIGP = IRHOP
   24 IF (IRHO-ISIG) 517,517,518
  518 ISTO = IRHO
      IRHO = ISIG
      ISIG = ISTO
  517 IF (IRHOP-ISIGP) 524,524,528
  528 ISTO = IRHOP
      IRHOP = ISIGP
      ISIGP = ISTO
  524 IF(IBUG2.GT.0) WRITE(IWRITE,61) IRHO,ISIG,IRHOP,ISIGP
      IF (NOVLPS .EQ. 0) GO TO 70
      MGO=MATCH(IA,IB,IC,ID)
      IF(MGO.EQ.2) GO TO 70
      IF(MGO.EQ.0) RETURN
      IF(IBACK3.EQ.0) RETURN
      GO TO 103
*
* --- CALCULATE COEFFICIENTS OF SLATER INTEGRALS
*
   70 CALL REDUCE(IRHO,ISIG,IRHOP,ISIGP,LESSEN)
      CALL ALLADD(IHSH,M3,M4,M5,M6,M7,M8,M9,M10,
     :M11,M12,M13,M14,M15,M16,M17,M18)
   75 CALL FANO(IRHO,ISIG,IRHOP,ISIGP)
      IF(LESSEN.NE.0) CALL MEREST(IRHO,ISIG,IRHOP,ISIGP)
      CALL PRNTWT(IRHO,ISIG,IRHOP,ISIGP)
      IF(NOVSTO.EQ.0) RETURN
  103 MU=MUSTO
      NU=NUSTO
      MUP=MUPSTO
      NUP=NUPSTO
      NOVLPS=NOVSTO
      LMU=LMUSTO
      LNU=LNUSTO
      LMUP=LMUPST
      LNUP=LNUPST
      IF(IBACK3.NE.0) GO TO 102
      RETURN
*
* === ONE INTERACTING SHELL SPECIFIED ON EACH SIDE. SUMMATION OVER OTHER
*
   20 IASTO=IA
      IBSTO=IB
      ICSTO=IC
      IDSTO=ID
      IRSTO=IRHO
      IRPSTO=IRHOP
      DO 125 K1=1,IHSH
      IF(NOVLPS.EQ.0) GO TO 39
      IA=IASTO
      IB=IBSTO
      IC=ICSTO
      ID=IDSTO
      IBACK1=0
      IBACK2=0
*
* --- STORE CURRENT DATA ON NON-ORTHOGONAL SUBSHELLS
*
      MUSTO=MU
      NUSTO=NU
      MUPSTO=MUP
      NUPSTO=NUP
      NOVSTO=NOVLPS
      LMUSTO=LMU
      LNUSTO=LNU
      LMUPST=LMUP
      LNUPST=LNUP
      IF(K1.EQ.MU) GO TO 40
      IF(K1.EQ.NU) GO TO 41
*
*    :SIG IS AN ORTHOGONAL SUBSHELL
*
   39 IF(NOSH1(K1)) 26,125,26
   26 ISIG=K1
      IF(NOSH2(K1)) 29,125,29
   29 ISIGP = K1
      IF(NOVLPS.EQ.0) GO TO 42
      MGO = MATCH(IA,IB,IC,ID)
      IF (MGO .EQ. 0) RETURN
      IF ( MGO .EQ. 1) GO TO 48
      GO TO 42
   40 IBACK1=1
      IF (NOVLPS .NE. 2 .OR. MUP .EQ. NUP) IBACK1 = 0
      ISIG = MU
      ISIGP = MUP
      IA=IA+1
      IC=IC+1
      GO TO 63
   46 IBACK1=0
      IF(IBACK3.EQ.1.AND.IRHO.NE.MU) GO TO 125
      ISIG=MU
      ISIGP=NUP
      IC=IC-1
      ID=ID+1
      GO TO 63
   41 IBACK2=1
      ISIG=NU
      ISIGP=NUP
      IB=IB+1
      ID=ID+1
      IF(IBACK3.EQ.1) THEN
         IF(IRHO.EQ.MU) GO TO 125
         GO TO 62
      ENDIF
      IF (MUP .EQ. NUP) GO TO 62
      GO TO 63
   62 IBACK2 = 0
      ISIG = NU
      ISIGP = MUP
      ID=ID-1
      IC=IC+1
   63 MGO = MATCH(IA,IB,IC,ID)
      IF (MGO .EQ. 0) RETURN
      IF (MGO .EQ. 2) GO TO 42
      GO TO 48
   47 IBACK2=0
      ISIG=NU
      ISIGP=MUP
      ID=ID-1
      IC=IC+1
      GO TO 63
   42 IRHO=IRSTO
      IRHOP=IRPSTO
*
*     ORTHOGONALITY OF THE ORBITALS REQUIRES THAT THE VARYING INTER-
*     ACTING SHELL BE THE SAME ON BOTH SIDES OF THE MATRIX ELEMENT
*
* --- IRHO.LE.ISIG,   IRHOP.LE.ISIGP
*
      IF(IRHO-ISIG) 27,127,227
  227 ISTO=IRHO
      IRHO=ISIG
      ISIG = ISTO
      GO TO 27
  127 IF(NOSH1(ISIG).LE.1) GO TO 48
   27 IF(IRHOP-ISIGP) 30,130,31
   31 ISTO=IRHOP
      IRHOP = ISIGP
      ISIGP = ISTO
      GO TO 30
  130 IF(NOSH2(ISIGP)-1) 48,48,30
   30 IF(IBUG2.GT.0) WRITE(IWRITE,61) IRHO,ISIG,IRHOP,ISIGP
      IF (K1 .EQ. IRSTO .OR. K1 .EQ. IRPSTO ) GO TO 71
      LC = 4*LJ(K1) + 2
      IF (NOSH1(K1) .NE. LC .OR. NOSH2(K1) .NE. LC) GO TO 71
      CALL CLSHEL(IRSTO,IRPSTO,K1)
      GO TO 48
*
* --- CALCULATE COEFFICIENTS OF SLATER INTEGRALS
*
   71 CALL REDUCE(IRHO,ISIG,IRHOP,ISIGP,LESSEN)
      CALL ALLADD(IHSH,M3,M4,M5,M6,M7,M8,M9,M10,
     :M11,M12,M13,M14,M15,M16,M17,M18)
   76 CALL FANO(IRHO,ISIG,IRHOP,ISIGP)
      IF(LESSEN.NE.0) CALL MEREST(IRHO,ISIG,IRHOP,ISIGP)
      CALL PRNTWT(IRHO,ISIG,IRHOP,ISIGP)
   48 IF(NOVSTO.EQ.0) GO TO 125
      MU=MUSTO
      NU=NUSTO
      MUP=MUPSTO
      NUP=NUPSTO
      NOVLPS=NOVSTO
      LMU=LMUSTO
      LNU=LNUSTO
      LMUP=LMUPST
      LNUP=LNUPST
      IRHO = IRSTO
      IRHOP = IRPSTO
      IF(IBACK1.EQ.1) GO TO 46
      IF(IBACK2.EQ.1) GO TO 47
  125 CONTINUE
      IF (NCLOSD.EQ.0) GO TO 128
      IF (NOVLPS.NE.0) THEN
         IA = IASTO
         IB = IBSTO
         IC = ICSTO
         ID = IDSTO
         MGO = MATCH(IA,IB,IC,ID)
         IF (MGO .NE. 2) GO TO 129
      END IF
      CALL CLSHEL(IRSTO,IRPSTO,IZERO)
129   IF (NOVSTO.NE.0) THEN
         MU = MUSTO
         NU = NUSTO
         MUP = MUPSTO
         NUP = NUPSTO
         NOVLPS = NOVSTO
         LMU = LMUSTO
         LNU = LNUSTO
         LMUP = LMUPST
         LNUP = LNUPST
         IRHO = IRSTO
         IRHOP = IRPSTO
      ENDIF
128   IF (IBACK3.EQ.1) GO TO 102
      RETURN
*
* === NO INTERACTING SHELLS SPECIFIED
*     SUMMATION OVER ALL POSSIBLE COMBINATIONS
*     IN THIS CASE, ORTHOGONALITY OF ORBITALS PRECLUDES ALL CASES
*     EXCEPT  IRHO=IRHOP    AND    ISIG=ISIGP
*
   19 DO 32 K1=1,IHSH
      IF(NOSH1(K1)) 36,32,36
   36 ISIG=K1
      DO 33 K2=1,K1
      IF(NOSH1(K2)) 37,33,37
   37 IRHO=K2
      IF(IRHO-ISIG) 50,51,50
   51 IF(NOSH1(ISIG)-1) 33,33,50
   50 IF(NOVLPS.EQ.0) GO TO 53
      MUSTO=MU
      NUSTO=NU
      MUPSTO=MUP
      NUPSTO=NUP
      NOVSTO=NOVLPS
      LMUSTO=LMU
      LNUSTO=LNU
      LMUPST=LMUP
      LNUPST=LNUP
      IF (IRHO.EQ.MU .AND. ISIG.EQ.MU) GO TO 81
      IF (IRHO.EQ.NU .AND. ISIG.EQ.NU) GO TO 82
      IF (IRHO.EQ.MU .AND. ISIG.EQ.NU) GO TO 83
      IF (IRHO.EQ.MU .OR. IRHO.EQ.NU) GO TO 84
      IF (ISIG.EQ.MU .OR. ISIG.EQ.NU) GO TO 85
      IRHOP = IRHO
      ISIGP = ISIG
      MGO = MATCH(IZERO,IZERO,IZERO,IZERO)
      IBACK = 3
   87 IF (MGO .EQ. 0) RETURN
      IF (MGO .EQ. 1) GO TO 69
      GO TO 58
   81 IBACK = 1
      IRHOP = MUP
      ISIGP = MUP
      MGO = MATCH(ITWO,IZERO,ITWO,IZERO)
      GO TO 87
   88 IBACK = 2
      IF (MUP .EQ. NUP) GO TO 33
      IRHOP = NUP
      ISIGP = NUP
      MGO = MATCH(ITWO,IZERO,IZERO,ITWO)
      GO TO 87
   89 IBACK = 3
      IRHOP = MUP
      ISIGP = NUP
      MGO = MATCH(ITWO,IZERO,IONE,IONE)
      GO TO 87
   82 IBACK = 4
      IRHOP = MUP
      ISIGP = MUP
      MGO = MATCH(IZERO,ITWO,ITWO,IZERO)
      GO TO 87
   90 IBACK = 5
      IF (MUP .EQ. NUP) GO TO 33
      IRHOP = NUP
      ISIGP = NUP
      MGO = MATCH(IZERO,ITWO,IZERO,ITWO)
      GO TO 87
   91 IBACK = 3
      IRHOP = MUP
      ISIGP = NUP
      MGO = MATCH(IZERO,ITWO,IONE,IONE)
      GO TO 87
   83 IBACK = 6
      IRHOP = MUP
      ISIGP = MUP
      MGO = MATCH(IONE,IONE,ITWO,IZERO)
      GO TO 87
   92 IBACK = 7
      IF (MUP .EQ. NUP) GO TO 33
      IRHOP = NUP
      ISIGP = NUP
      MGO = MATCH(IONE,IONE,IZERO,ITWO)
      GO TO 87
   93 IBACK = 3
      IRHOP = MUP
      ISIGP = NUP
      MGO = MATCH(IONE,IONE,IONE,IONE)
      GO TO 87
   84 ISIGP = ISIG
      IF (IRHO .EQ. NU .AND. MU.NE.NU) GO TO 94
      IBACK = 8
      IRHOP = MUP
      MGO = MATCH(IONE,IZERO,IONE,IZERO)
   96 IF (IRHOP .LT. ISIGP) GO TO 87
      IST = ISIGP
      ISIGP = IRHOP
      IRHOP = IST
      GO TO 87
   95 ISIGP = ISIG
      IBACK = 3
      IF (MUP .EQ. NUP) GO TO 33
      IRHOP = NUP
      MGO = MATCH(IONE,IZERO,IZERO,IONE)
      GO TO 96
   94 IBACK = 9
      IRHOP = MUP
      MGO = MATCH(IZERO,IONE,IONE,IZERO)
      GO TO 96
   97 ISIGP = ISIG
      IF (MUP .EQ. NUP) GO TO 33
      IBACK = 3
      IRHOP = NUP
      MGO = MATCH(IZERO,IONE,IZERO,IONE)
      GO TO 96
   85 IRHOP = IRHO
      IF (ISIG .EQ. NU .AND. MU .NE. NU) GO TO 98
      IBACK = 10
      ISIGP = MUP
      MGO = MATCH(IONE,IZERO,IONE,IZERO)
      GO TO 96
   99 IRHOP = IRHO
      IF (MUP .EQ. NUP) GO TO 33
      IBACK = 3
      ISIGP = NUP
      MGO = MATCH(IONE,IZERO,IZERO,IONE)
      GO TO 96
   98 IBACK = 11
      ISIGP = MUP
      MGO = MATCH(IZERO,IONE,IONE,IZERO)
      GO TO 96
   79 IRHOP = IRHO
      IF (MU .EQ. NUP) GO TO 33
      IBACK = 3
      ISIGP = NUP
      MGO = MATCH(IZERO,IONE,IZERO,IONE)
      GO TO 96
   53 IRHOP=IRHO
      ISIGP=ISIG
   58 IF(IBUG2.GT.1) WRITE(IWRITE,61) IRHO,ISIG,IRHOP,ISIGP
*
* --- CALCULATE COEFFICIENTS OF SLATER INTEGRALS
*
   74 IF(ICHOP(K1).EQ.1.OR.ICHOP(K2).EQ.1) GO TO 34
      CALL REDUCE(IRHO,ISIG,IRHOP,ISIGP,LESSEN)
      CALL ALLADD(IHSH,M3,M4,M5,M6,M7,M8,M9,M10,
     :M11,M12,M13,M14,M15,M16,M17,M18)
   77 CALL FANO(IRHO,ISIG,IRHOP,ISIGP)
      IF(LESSEN.NE.0) CALL MEREST(IRHO,ISIG,IRHOP,ISIGP)
      CALL PRNTWT(IRHO,ISIG,IRHOP,ISIGP)
      GO TO 69
   34 IF (IFULL.NE.0 .OR. IALL.NE.0) CALL USEEAV(IRHO,ISIG)
   69 IF(NOVSTO.EQ. 0) GO TO 33
      MU=MUSTO
      NU=NUSTO
      MUP=MUPSTO
      NUP=NUPSTO
      NOVLPS=NOVSTO
      LMU=LMUSTO
      LNU=LNUSTO
      LMUP=LMUPST
      LNUP=LNUPST
      GO TO (88,89,33,90,91,92,93,95,97,99,79) IBACK
   33 CONTINUE
   32 CONTINUE
      IF (NCLOSD .EQ. 0 ) RETURN
      ISIG = 0
      CALL CLSHEL(ISIG,ISIG,IZERO)
    1 RETURN
*
* --- SET CONSTANTS USEFUL IN INNER SUBROUTINES
*
      END
*
*     ------------------------------------------------------------------
*       S C R A P
*     ------------------------------------------------------------------
*
      SUBROUTINE SCRAP(IX,IRS,IR1)
      PARAMETER(KFL2=12)
      DIMENSION IX(KFL2,3)
*
* === Deletes row IR1 from the coupling array I1
*
      IF(IR1.EQ.IRS) GO TO 3
      IR2=IR1+1
      DO 1 I=IR2,IRS
      I1=I-1
      DO 2 K=1,3
      IX(I1,K)=IX(I,K)
    2 CONTINUE
    1 CONTINUE
      GO TO 4
    3 IX(IRS-1,3)=IX(IRS,3)
    4 IRS=IRS-1
      RETURN
      END
*
*     ------------------------------------------------------------------
*       S E T J 1
*     ------------------------------------------------------------------
*
      SUBROUTINE SETJ1(K,IRHO,ISIG,IRHOP,ISIGP,ITST1,ITST2,ITST3,ITST4,
     :                 K1,K2,K3,K4)
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      PARAMETER (NWD=30,NCD=100,NORD=NWD*(NWD-1)/2)
*
      COMMON/INFORM/IREAD,IWRITE,IOUT,ISC0,ISC1,ISC2,ISC3,JSC(3),IALL
      COMMON/DEBUG/IBUG1,IBUG2,IBUG3,NBUG6,NBUG7,IFULL
      COMMON/MEDEFN/IHSH,NJ(10),LJ(10),NOSH1(10),NOSH2(10),J1QN1(19,3),
     :     J1QN2(19,3),IJFUL(10)
      COMMON/MVALUE/M1,M2,M3,M4,M5,M6,M7,M8,M9,M10,M11,M12,M13,M14,M15,
     :     M16,M17,M18,M19,M20
      COMMON/NJLJ/NRHO,LRHO,NSIG,LSIG,NRHOP,LRHOP,NSIGP,LSIGP
      COMMON/INTERM/J1BAR1(10,3),J1BAR2(10,3),J1TLD1(3),J1TLD2(3)
      PARAMETER(KFL1=60,KFL2=12)
      LOGICAL FREE
      COMMON/COUPLE/NJ1S,NJ23S,J1(KFL1),J2(KFL2,3),J3(KFL2,3),FREE(KFL1)
      COMMON/OVRLAP/MU,NU,MUP,NUP,NONORT,NOVLPS,IROWMU,IROWNU,ICOLMU,
     : ICOLNU,NORTH,IORDER,NCALLS,LMU,LNU,LMUP,LNUP,JMU,JNU,JMUP,JNUP,
     :     IORTH(NORD)
*
* === SETS J1 AND FREE ARRAYS FOR DIRECT INTEGRAL CALLS OF NJGRAF
*
   14 FORMAT('  J1:  ',36I3)
   15 FORMAT('FREE:  ',36L3)
      DO 1 J=1,IHSH
   17 J1(J)=J1BAR1(J,K)
      IF (J.EQ.MUP .OR. J.EQ.NUP) J1(J) = J1BAR2(J,K)
    1 CONTINUE
      DO 2 J=M4,M6
      J1(J)=J1QN1(J,K)
    2 CONTINUE
      DO 3 J=M7,M8
      J1(J)=J1QN2(J-M3,K)
    3 CONTINUE
      J1(M10)=J1QN1(ISIG,K)
      J1(M12)=J1QN2(ISIGP,K)
      IF(M1) 4,5,4
    4 J1(M9)=J1QN1(IRHO,K)
      GO TO 6
    5 J1(M9)=J1TLD1(K)
    6 IF(M2) 7,8,7
    7 J1(M11)=J1QN2(IRHOP,K)
      GO TO 9
    8 J1(M11)=J1TLD2(K)
*
*     K=2  IMPLIES ANGULAR TERMS   ,   K=3  IMPLIES SPIN TERMS
*
    9 IF(K-2) 11,11,10
   10 J1(M13)=2
      J1(M14)=2
      MLIMIT=M14
      NJ1S=M14
      NJ23S=M5
      GO TO 12
   11 J1(M13)=2*LRHO+1
      J1(M14)=2*LSIG+1
      J1(M15)=2*LRHOP+1
      J1(M16)=2*LSIGP+1
      MLIMIT=M16
      NJ1S=M17
      NJ23S=NJ23S+1
*
*     PRINT-OUT OF VALUES IN NJSYM  IF IBUG3=1
*
   12 IF(IBUG2.GT.0 .AND. IBUG3.NE.1)
     :  WRITE(IWRITE,14)(J1(J),J=1,NJ1S)
*
*     If ITST1.eq.0 or ITST2.eq.0 or ITST3.eq.0 or ITST4.eq.0 then
*     NJGRAF is being called, so set the elements of the FREE array.
*
      IF(ITST1.EQ.0.OR.ITST2.EQ.0.OR.ITST3.EQ.0.OR.ITST4.EQ.0) THEN
          DO 20 J=1,MLIMIT
              FREE(J)=.FALSE.
   20     CONTINUE
*
          IF((K1.NE.1).AND.(IRHO.NE.MUP).AND.(IRHO.NE.NUP))
     :    FREE(IRHO)=.TRUE.
          IF((K3.NE.1).AND.(IRHOP.EQ.MUP.OR.IRHOP.EQ.NUP))
     :    FREE(IRHOP)=.TRUE.
          IF(M1.EQ.0) THEN
              IF(K2.NE.1) FREE(M9)=.TRUE.
          ELSE
              IF((K2.NE.1).AND.(ISIG.NE.MUP).AND.(ISIG.NE.NUP))
     :        FREE(ISIG)=.TRUE.
          ENDIF
          IF(M2.EQ.0) THEN
              IF(K4.NE.1) FREE(M11)=.TRUE.
          ELSE
              IF((K4.NE.1).AND.(ISIGP.EQ.MUP.OR.ISIGP.EQ.NUP))
     :        FREE(ISIGP)=.TRUE.
          ENDIF
*
*     Print-out of values in NJGRAF if IBUG2>0.
*
          IF(IBUG2.GT.0.AND.IBUG3.NE.1)
     :      WRITE(IWRITE,15) (FREE(J),J=1,MLIMIT)
*
      ENDIF
*
      END
*
*     ------------------------------------------------------------------
*       S E T U P
*     ------------------------------------------------------------------
*
      SUBROUTINE SETUP(JA,JB)
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      PARAMETER (NWD=30,NCD=100,NORD=NWD*(NWD-1)/2)
*
      COMMON/INFORM/IREAD,IWRITE,IOUT,ISC0,ISC1,ISC2,ISC3,JSC(3),IALL
      COMMON/DEBUG/IBUG1,IBUG2,IBUG3,NBUG6,NBUG7,IFULL
      COMMON/STATES/NCFG,MAXORB,IAJCMP(NWD),LJCOMP(NWD),
     :NJCOMP(NWD),NOCCSH(NCD),NELCSH(5,NCD),NOCORB(5,NCD),J1QNRD(9,NCD)
      COMMON/MEDEFN/IHSH,NJ(10),LJ(10),NOSH(10,2),J1QN(19,3,2),IJFUL(10)
      COMMON/OVRLAP/MU,NU,MUP,NUP,NONORT,NOVLPS,IROWMU,IROWNU,ICOLMU,
     : ICOLNU,NORTH,IORDER,NCALLS,LMU,LNU,LMUP,LNUP,JMU,JNU,JMUP,JNUP,
     :     IORTH(NORD)
*
*     NOTICE THE DIFFERENT NAMES IN THE COMMON BLOCK MEDEFN  -  WE
*      STORE NOSH1(I=1,10) AS NOSH((I=1,10),1) AND NOSH2(I=1,10) AS
*     NOSH((I=1,10),2)   AND USE THE FACT THAT NOSH1 AND NOSH2 WILL THEN
*     BE EQUIVALENT TO THE SINGLE 2-DIMENSIONAL ARRAY NOSH.  SIMILARLY
*     FOR J1QN
*
* === GENERATES THE ARRAYS  NJ,LJ - DEFINING THE QUANTUM NUMBERS OF THE
*     SHELLS,   NOSH - DEFINING THE OCCUPATION OF THE SHELLS,  J1QN -
*     DEFINING THE COUPLING OF THE SHELLS,   FOR EACH OF THE TWO
*     CONFIGURATIONS CONSIDERED.    ONLY THOSE SHELLS OCCURRING IN AT
*     LEAST ONE CONFIGURATION ARE INCLUDED.
*                   AT LEAST TWO SHELLS MUST BE CONSIDERED OCCUPIED.
*     THUS (1S)**2    HELIUM  MUST BE TREATED AS ,E.G., (1S)**2(2S)**0
*     THE SIZE OF THE ARRAYS HERE CALCULATED IS ARRANGED TO BE NO
*     GREATER THAN IS NECESSARY TO INCLUDE ALL ORBITALS WHICH ARE
*     DEEMED TO BE OCCUPIED IN EITHER OR BOTH OF THE CONFIGURATIONS
*     JA,JB
*
* --- INITIALIZE BASIC QUANTITIES - (I1+1) RUNS OVER 1,MAXORB,  IHSH IS
*     THE CURRENT VALUE OF THE HIGHEST OCCUPIED SHELL YET CONSIDERED,
*     WHILE I2HSH=2*IHSH-1
*
      MU=0
      NU=0
      MUP=0
      NUP=0
      I1=0
      IHSH=0
      I2HSH=-1
      IA=NOCCSH(JA)
      IB=NOCCSH(JB)
*
* --- TEST ON WHETHER LIMIT OF I1 HAS BEEN REACHED
*
    1 IF(I1-MAXORB) 101,100,100
*
* --- INCREASE BASIC QUANTITIES
*
  101 I1=I1+1
      I3=IHSH+1
      I5=I2HSH+I3
*
* --- IS THE SHELL I1 OCCUPIED IN JA
*
      DO 2 J=1,IA
      IF(I1-NOCORB(J,JA)) 2,3,2
    2 CONTINUE
      NA=1
      GO TO 4
    3 NA=2
      J1=J
*
* --- IS THE SHELL I1 OCCUPIED IN JB
*
    4 DO 5 J=1,IB
      IF(I1-NOCORB(J,JB)) 5,6,5
    5 CONTINUE
      NB=1
      GO TO 7
    6 NB=2
      J2=J
*
*     IF THE SHELL I1 IS NOT OCCUPIED IN EITHER JA OR JB, IGNORE THE
*     SHELL, DO NOT INCREASE IHSH, AND CONSIDER NEXT SHELL BY INCREASING
*     I1
*
    7 IF(NA-1) 8,8,9
    8 IF(NB-1) 1,1,9
*
* --- IF THE SHELL I1 IS OCCUPIED IN EITHER JA OR JB -
*     (1)   IF IHSH.GT.1, THEN ALREADY AT LEAST TWO SHELLS AND THE
*     RESULTING COUPLINGS HAVE BEEN STORED. WE MUST THUS MAKE ROOM FOR
*     THE QUANTUM NUMBERS OF THIS NEW SHELL BETWEEN THE QUANTUM NUMBERS
*     OF THE PREVIOUS SHELLS AND THE QUANTUM NUMBERS OF THE INTERMEDIATE
*     COUPLINGS OF THE CONFIGURATIONS.  THUS THE LATTER SET ARE =MOVED
*     ALONG= TO MAKE ROOM FOR THE NEW SHELL
*     (2)   IF IHSH.LE.1, THERE ARE NO INTERMEDIATE COUPLING QUANTUM
*     NUMBERS, AND SO THERE IS NOTHING TO MOVE
*
    9 IF(IHSH-1) 11,11,10
   10 DO 12 I=1,2
      DO 13 J=I3,I2HSH
      I4=I5-J
      DO 14 K=1,3
      J1QN(I4+1,K,I)=J1QN(I4,K,I)
   14 CONTINUE
   13 CONTINUE
   12 CONTINUE
   11 IHSH=I3
      I2HSH=I2HSH+2
      NC=NA
      I=1
      IC=J1
      JC=JA
*
* --- FIRST CONSIDER THE L.H.S. (I=1) OF THE MATRIX ELEMENT.  NC=1 MEANS
*     UNOCCUPIED, REPRESENTED BY A DUMMY SINGLET S SHELL, AND THE
*     ADDITIONAL SET OF COUPLING QUANTUM NUMBERS WILL BE THE SAME AS THE
*     LAST SET OF COUPLING QUANTUM NUMBERS ALREADY OBTAINED.
*     NC=2 MEANS OCCUPIED.  THEN ALL THE NEW QUANTUM NUMBERS (BOTH FOR
*     THE SHELL AND FOR THE COUPLING OF THIS SHELL TO THE RESULTANT OF
*     THE PREVIOUS ONES) ARE DEFINED IN THE CORRESPONDING J1QNRD ARRAY.
*     NOSH - THE NUMBER OF ELECTRONS IN THIS SHELL, IS DEFINED BY THE
*     APPROPRIATE ENTRY IN NELCSH .  THE R.H.S. IS THEN CONSIDERED
*     SIMILARLY (I=2)
*
   25 GO TO (15,16),NC
   15 NOSH(IHSH,I)=0
      J1QN(IHSH,1,I)=0
      J1QN(IHSH,2,I)=1
      J1QN(IHSH,3,I)=1
      IF(IHSH-2) 22,18,19
   18 J1QN(3,1,I)=0
      J1QN(3,2,I)=J1QN(1,2,I)
      J1QN(3,3,I)=J1QN(1,3,I)
      GO TO 22
   19 DO 27 K=1,3
      J1QN(I2HSH,K,I)=J1QN(I2HSH-1,K,I)
   27 CONTINUE
      GO TO 22
   16 NOSH(IHSH,I)=NELCSH(IC,JC)
      IF(NOVLPS.EQ.0) GO TO 33
      GO TO (31,32),I
   31 IF(IC.EQ.JMU) MU=IHSH
      IF(IC.EQ.JNU) NU=IHSH
      GO TO 33
   32 IF(IC.EQ.JMUP) MUP=IHSH
      IF(IC.EQ.JNUP) NUP=IHSH
   33 JD = J1QNRD(IC,JC)
      J1QN(IHSH,1,I)=MOD(JD,64)
      JD = JD/64
      J1QN(IHSH,2,I) = MOD(JD,64)
      J1QN(IHSH,3,I) = JD/64
*
*     IS THIS THE FIRST OCCUPIED SHELL OF EITHER CONFIGURATION. IF SO,
*     THEN THERE ARE NO INTERMEDIATE COUPLINGS TO CONSIDER AT THIS STAGE
*
      IF(IHSH .GT. 1) THEN
*
*     IS THIS THE FIRST OCCUPIED SHELL OF THIS CONFIGURATION, THOUGH NOT
*     THE FIRST OF THE OTHER CONFIGURATION.  IF SO, THE INTERMEDIATE
*     COUPLING FORMED HAS THE SAME  L,S  VALUES AS THIS OCCUPIED SHELL,
*     SINCE WE COUPLE THE SHELL TO A DUMMY SINGLET S.
*
         IF(IC .LE.1) THEN 
            I2 = 1
         ELSE
            I2 = NOCCSH(JC)+IC-1
         END IF
         JD = J1QNRD(I2,JC)
         IF (IC .LE. 1) THEN
            J1QN(I2HSH,1,I) = 0
         ELSE
            J1QN(I2HSH,1,I) = MOD(JD,64)
         END IF
         JD = JD/64
         J1QN(I2HSH,2,I) = MOD(JD, 64)
         J1QN(I2HSH,3,I) = JD/64
      END IF
*
*     SENIORITY SET (ARBITRARILY) ZERO FOR INTERMEDIATE COUPLING
*
   22 IF(I-2) 23,24,24
   23 NC=NB
      I=2
      IC=J2
      JC=JB
      GO TO 25
*
* --- SET THE NJ AND LJ VALUES OF THE OCCUPIED SHELLS
*
   24 NJ(IHSH)=NJCOMP(I1)
      IJFUL(IHSH)=I1
      LJ(IHSH)=LJCOMP(I1)
*
* --- RETURN TO 1  TO SEE IF MAXORB HAS BEEN REACHED
*
      GO TO 1
  100 CONTINUE
      END
*
*     ------------------------------------------------------------------
*       U S E E A V
*     ------------------------------------------------------------------
*
      SUBROUTINE USEEAV(IRHO,ISIG)
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      PARAMETER (NWD=30,NCD=100,NORD=NWD*(NWD-1)/2)
*
      COMMON/INFORM/IREAD,IWRITE,IOUT,ISC0,ISC1,ISC2,ISC3,JSC(3),IALL
      COMMON/DEBUG/IBUG1,IBUG2,IBUG3,NBUG6,NBUG7,IFULL
      COMMON/STATES/NCFG,MAXORB,IAJCMP(NWD),LJCOMP(NWD),
     :NJCOMP(NWD),NOCCSH(NCD),NELCSH(5,NCD),NOCORB(5,NCD),J1QNRD(9,NCD)
      COMMON/ENAV/NINTS,KVALUE(15),COEFCT(15)
      COMMON/MEDEFN/IHSH,NJ(10),LJ(10),NOSH1(10),NOSH2(10),J1QN1(19,3),
     :     J1QN2(19,3),IJFUL(10)
      COMMON/DIAGNL/IDIAG,JA,JB
*
*     DETERMINE THE INTERACTION ENERGY
*
      JRHO=IJFUL(IRHO)
      JSIG=IJFUL(ISIG)
      N1=NOSH1(IRHO)
      N2=NOSH2(ISIG)
      M1=ISIG-IRHO
      IZERO=0
      ZERO=0.D0
      IF(M1.EQ.0) GO TO 1
      IEQUIV=2
      AC2=DFLOAT(N1*N2)
      GO TO 2
    1 IEQUIV=1
      AC2=DFLOAT(N1*(N1-1)/2)
    2 LA=LJ(IRHO)
      LB=LJ(ISIG)
      CALL INTACT(LA,LB,IEQUIV)
      IF(IBUG2-1) 4,7,7
*
*     PRINT OUT RESULTS AS IN SUBROUTINE PRNTWT
*
    7  WRITE(IWRITE,10) IAJCMP(JRHO),IAJCMP(JSIG),
     : IAJCMP(JRHO),IAJCMP(JSIG)
    4 IF (IFULL .NE. 0)
     :  WRITE(IWRITE,6)AC2,AC2,ZERO,IZERO,IAJCMP(JRHO),IAJCMP(JSIG)
      IF (IALL .NE. 0) CALL SAVE(1,AC2,IZERO,0,JRHO,0,JSIG,JA,JB,0)
      IF(NINTS.EQ.0) RETURN
      DO 8 N=1,NINTS
      ZA=AC2*COEFCT(N)
      K=KVALUE(N)
      IF(IEQUIV.EQ.1) GO TO 9
      IF (IFULL .NE. 0)
     :   WRITE(IWRITE,11) ZA,ZA,ZERO,K,IAJCMP(JRHO),IAJCMP(JSIG)
      IF (IALL .NE. 0) CALL SAVE(2,ZA,K,0,JRHO,0,JSIG,JA,JB,0)
      GO TO 8
    9 IF (IFULL .NE. 0)
     :   WRITE(IWRITE,6) ZA,ZA,ZERO,K,IAJCMP(JRHO),IAJCMP(JSIG)
      IF (IALL .NE. 0) CALL SAVE(1,ZA,K,0,JRHO,0,JSIG,JA,JB,0)
    8 CONTINUE
    6 FORMAT(3F14.6,4X,1HF,I2,1H(,A3,1H,,A3,1H))
   10 FORMAT(//23H INTERACTING SHELLS ARE,6X,6H RHO =,A3,6X,6H SIG =,A3,
     :6X,7H RHOP =,A3,6X,7H SIGP =,A3//)
   11 FORMAT(3F14.6,4X,1HG,I2,1H(,A3,1H,,A3,1H))
  100 FORMAT(F14.8,1HF,I2,1H(,A3,I3,1H,,A3,I3,1H))
  101 FORMAT(F14.8,1HG,I2,1H(,A3,I3,1H,,A3,I3,1H))
      RETURN
      END
